mirror of
https://github.com/TREX-CoE/qmckl.git
synced 2024-08-25 06:31:45 +02:00
Merge branch 'master' into add-python-api
This commit is contained in:
commit
9711f2699c
116
.github/workflows/test-build.yml
vendored
116
.github/workflows/test-build.yml
vendored
@ -2,9 +2,7 @@ name: test-build
|
||||
|
||||
on:
|
||||
push:
|
||||
branches: [ master ]
|
||||
pull_request:
|
||||
branches: [ master ]
|
||||
|
||||
jobs:
|
||||
x86_ubuntu:
|
||||
@ -78,60 +76,60 @@ jobs:
|
||||
run: make python-test
|
||||
working-directory: _build
|
||||
|
||||
x86_macos:
|
||||
|
||||
runs-on: macos-latest
|
||||
name: x86 MacOS latest
|
||||
|
||||
steps:
|
||||
- uses: actions/checkout@v2
|
||||
- name: install dependencies
|
||||
run: brew install emacs hdf5 automake pkg-config
|
||||
|
||||
- name: Symlink gfortran (macOS)
|
||||
if: runner.os == 'macOS'
|
||||
run: |
|
||||
# make sure gfortran is available
|
||||
# https://github.com/actions/virtual-environments/issues/2524
|
||||
# https://github.com/cbg-ethz/dce/blob/master/.github/workflows/pkgdown.yaml
|
||||
sudo ln -s /usr/local/bin/gfortran-10 /usr/local/bin/gfortran
|
||||
sudo mkdir /usr/local/gfortran
|
||||
sudo ln -s /usr/local/Cellar/gcc@10/*/lib/gcc/10 /usr/local/gfortran/lib
|
||||
gfortran --version
|
||||
|
||||
- name: Install the latest TREXIO from the GitHub clone
|
||||
run: |
|
||||
git clone https://github.com/TREX-CoE/trexio.git
|
||||
cd trexio
|
||||
./autogen.sh
|
||||
./configure --prefix=${PWD}/_install --enable-silent-rules
|
||||
make -j 4
|
||||
make install
|
||||
|
||||
- name: Test TREXIO
|
||||
run: make -j 4 check
|
||||
working-directory: trexio
|
||||
|
||||
- name: Archive TREXIO test log file
|
||||
if: failure()
|
||||
uses: actions/upload-artifact@v2
|
||||
with:
|
||||
name: test-report-trexio-macos
|
||||
path: trexio/test-suite.log
|
||||
|
||||
- name: Build QMCkl
|
||||
run: |
|
||||
export PKG_CONFIG_PATH=${PWD}/trexio/_install/lib/pkgconfig:$PKG_CONFIG_PATH
|
||||
./autogen.sh
|
||||
./configure CC=gcc-10 FC=gfortran-10 --enable-silent-rules --enable-debug
|
||||
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-macos
|
||||
path: test-suite.log
|
||||
# x86_macos:
|
||||
#
|
||||
# runs-on: macos-latest
|
||||
# name: x86 MacOS latest
|
||||
#
|
||||
# steps:
|
||||
# - uses: actions/checkout@v2
|
||||
# - name: install dependencies
|
||||
# run: brew install emacs hdf5 automake pkg-config
|
||||
#
|
||||
# - name: Symlink gfortran (macOS)
|
||||
# if: runner.os == 'macOS'
|
||||
# run: |
|
||||
# # make sure gfortran is available
|
||||
# # https://github.com/actions/virtual-environments/issues/2524
|
||||
# # https://github.com/cbg-ethz/dce/blob/master/.github/workflows/pkgdown.yaml
|
||||
# sudo ln -s /usr/local/bin/gfortran-10 /usr/local/bin/gfortran
|
||||
# sudo mkdir /usr/local/gfortran
|
||||
# sudo ln -s /usr/local/Cellar/gcc@10/*/lib/gcc/10 /usr/local/gfortran/lib
|
||||
# gfortran --version
|
||||
#
|
||||
# - name: Install the latest TREXIO from the GitHub clone
|
||||
# run: |
|
||||
# git clone https://github.com/TREX-CoE/trexio.git
|
||||
# cd trexio
|
||||
# ./autogen.sh
|
||||
# ./configure --prefix=${PWD}/_install --enable-silent-rules
|
||||
# make -j 4
|
||||
# make install
|
||||
#
|
||||
# - name: Test TREXIO
|
||||
# run: make -j 4 check
|
||||
# working-directory: trexio
|
||||
#
|
||||
# - name: Archive TREXIO test log file
|
||||
# if: failure()
|
||||
# uses: actions/upload-artifact@v2
|
||||
# with:
|
||||
# name: test-report-trexio-macos
|
||||
# path: trexio/test-suite.log
|
||||
#
|
||||
# - name: Build QMCkl
|
||||
# run: |
|
||||
# export PKG_CONFIG_PATH=${PWD}/trexio/_install/lib/pkgconfig:$PKG_CONFIG_PATH
|
||||
# ./autogen.sh
|
||||
# ./configure CC=gcc-10 FC=gfortran-10 --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-macos
|
||||
# path: test-suite.log
|
||||
|
56
Makefile.am
56
Makefile.am
@ -93,6 +93,36 @@ html-local: $(htmlize_el) $(dist_html_DATA)
|
||||
text: $(htmlize_el) $(dist_text_DATA)
|
||||
doc: html text
|
||||
|
||||
setup_py = $(srcdir)/python/setup.py
|
||||
process_header_py = $(srcdir)/python/src/process_header.py
|
||||
test_py = $(srcdir)/python/test/test_api.py
|
||||
qmckl_i = $(srcdir)/python/src/qmckl.i
|
||||
numpy_i = $(srcdir)/python/src/numpy.i
|
||||
qmckl_wrap_c = python/src/qmckl_wrap.c
|
||||
qmckl_include_i = python/src/qmckl_include.i
|
||||
qmckl_py = python/src/qmckl.py
|
||||
|
||||
dist_python_DATA = $(setup_py) $(qmckl_py) $(qmckl_wrap_c) \
|
||||
$(srcdir)/python/pyproject.toml \
|
||||
$(srcdir)/python/requirements.txt \
|
||||
$(srcdir)/python/README.md
|
||||
|
||||
python-install: $(qmckl_h) $(lib_LTLIBRARIES) $(dist_python_DATA)
|
||||
$(MKDIR_P) python/src
|
||||
cd python ; \
|
||||
[[ ! -f pyproject.toml ]] && \
|
||||
cp $(abs_srcdir)/python/{pyproject.toml,requirements.txt,README.md,setup.py} . ; \
|
||||
cp src/qmckl.py . ; \
|
||||
export QMCKL_INCLUDEDIR="$(prefix)/include" ; \
|
||||
export QMCKL_LIBDIR="$(prefix)/lib" ; \
|
||||
pip install .
|
||||
|
||||
python-test: $(test_py)
|
||||
cd $(abs_srcdir)/python/test/ && \
|
||||
python test_api.py
|
||||
|
||||
.PHONY: python-test python-install cppcheck
|
||||
|
||||
|
||||
if QMCKL_DEVEL
|
||||
|
||||
@ -170,15 +200,6 @@ cppcheck.out: $(qmckl_h)
|
||||
--language=c --std=c99 -rp --platform=unix64 \
|
||||
-I$(srcdir)/include -I$(top_builddir)/include *.c *.h 2>../$@
|
||||
|
||||
setup_py = $(srcdir)/python/setup.py
|
||||
process_header_py = $(srcdir)/python/src/process_header.py
|
||||
test_py = $(srcdir)/python/test/test_api.py
|
||||
qmckl_i = $(srcdir)/python/src/qmckl.i
|
||||
numpy_i = $(srcdir)/python/src/numpy.i
|
||||
qmckl_wrap_c = python/src/qmckl_wrap.c
|
||||
qmckl_include_i = python/src/qmckl_include.i
|
||||
qmckl_py = python/qmckl/qmckl.py
|
||||
|
||||
$(qmckl_include_i): $(qmckl_h) $(process_header_py)
|
||||
$(MKDIR_P) python/src
|
||||
python $(process_header_py) $(qmckl_h)
|
||||
@ -186,28 +207,13 @@ $(qmckl_include_i): $(qmckl_h) $(process_header_py)
|
||||
|
||||
|
||||
$(qmckl_py): $(qmckl_i) $(qmckl_include_i)
|
||||
swig -Iinclude -Ipython/src -python -py3 -builtin -o $(qmckl_wrap_c) $(qmckl_i)
|
||||
$(srcdir)/tools/missing swig -Iinclude -Ipython/src -python -py3 -builtin -o $(qmckl_wrap_c) $(qmckl_i) || :
|
||||
|
||||
$(qmckl_wrap_c): $(qmckl_py)
|
||||
|
||||
python-install: $(qmckl_h) $(qmckl_i) $(setup_py) $(qmckl_py) $(qmckl_wrap_c)
|
||||
$(MKDIR_P) python/src
|
||||
cd python ; \
|
||||
[[ ! -f pyproject.toml ]] && \
|
||||
cp $(abs_srcdir)/python/{pyproject.toml,requirements.txt,README.md,setup.py} . ; \
|
||||
cp src/qmckl.py . ; \
|
||||
export QMCKL_INCLUDEDIR="$(prefix)/include" ; \
|
||||
export QMCKL_LIBDIR="$(prefix)/lib" ; \
|
||||
pip install .
|
||||
|
||||
python-test: $(test_py)
|
||||
cd $(abs_srcdir)/python/test/ && \
|
||||
python test_api.py
|
||||
|
||||
CLEANFILES += $(qmckl_wrap_c) \
|
||||
$(qmckl_include_i) \
|
||||
$(qmckl_py)
|
||||
|
||||
.PHONY: cppcheck python-test python-install
|
||||
|
||||
endif
|
||||
|
@ -1,5 +1,7 @@
|
||||
# QMCkl: Quantum Monte Carlo Kernel Library
|
||||
|
||||
<img src="https://trex-coe.eu/sites/default/files/styles/responsive_no_crop/public/2022-01/QMCkl%20code.png?itok=UvOUClA5" width=200>
|
||||
|
||||
![Build Status](https://github.com/TREX-CoE/qmckl/workflows/test-build/badge.svg?branch=master)
|
||||
|
||||
The domain of quantum chemistry needs a library in which the main
|
||||
|
@ -1,5 +1,5 @@
|
||||
#!/bin/bash
|
||||
|
||||
export srcdir="."
|
||||
python ${srcdir}/tools/build_makefile.py
|
||||
python3 ${srcdir}/tools/build_makefile.py
|
||||
autoreconf -i -Wall --no-recursive
|
||||
|
169
configure.ac
169
configure.ac
@ -35,7 +35,7 @@
|
||||
|
||||
AC_PREREQ([2.69])
|
||||
|
||||
AC_INIT([qmckl],[0.1.1],[https://github.com/TREX-CoE/qmckl/issues],[],[https://trex-coe.github.io/qmckl/index.html])
|
||||
AC_INIT([qmckl],[0.2.1],[https://github.com/TREX-CoE/qmckl/issues],[],[https://trex-coe.github.io/qmckl/index.html])
|
||||
AC_CONFIG_AUX_DIR([tools])
|
||||
AM_INIT_AUTOMAKE([subdir-objects color-tests parallel-tests silent-rules 1.11])
|
||||
|
||||
@ -93,6 +93,7 @@ AC_PROG_F77
|
||||
m4_version_prereq([2.70],[], [AC_PROG_CC_C99])
|
||||
AS_IF([test "$ac_cv_prog_cc_c99" = "no"], [AC_MSG_ERROR([The compiler does not support C99])])
|
||||
AC_PROG_CC_C_O
|
||||
AM_PROG_CC_C_O
|
||||
AC_PROG_FC
|
||||
AC_PROG_FC_C_O
|
||||
AC_FC_PP_DEFINE
|
||||
@ -212,6 +213,9 @@ esac
|
||||
|
||||
case $CC in
|
||||
|
||||
*gcc*)
|
||||
CFLAGS="$CFLAGS -fPIC"
|
||||
;;
|
||||
*nvc*)
|
||||
CFLAGS="$CFLAGS -fPIC"
|
||||
;;
|
||||
@ -224,6 +228,131 @@ AS_IF([test "$HAVE_HPC" = "yes"], [
|
||||
AC_DEFINE([HAVE_HPC], [1], [If defined, activate HPC routines])
|
||||
])
|
||||
|
||||
# Enable Verificarlo tests
|
||||
AC_ARG_ENABLE([vfc_ci],
|
||||
[ --enable-vfc_ci Build the library with vfc_ci support],
|
||||
[case "${enableval}" in
|
||||
yes) vfc_ci=true && FCFLAGS="-D VFC_CI $FCFLAGS" && CFLAGS="-D VFC_CI $CFLAGS";;
|
||||
no) vfc_ci=false ;;
|
||||
*) AC_MSG_ERROR([bad value ${enableval} for --enable_vfc_ci]) ;;
|
||||
esac],[vfc_ci=false])
|
||||
AM_CONDITIONAL([VFC_CI], [test x$vfc_ci = xtrue])
|
||||
|
||||
if test "$FC" = "verificarlo-f"; then
|
||||
AC_MSG_NOTICE(verificarlo-f detected)
|
||||
# Arguments order is important here
|
||||
FCFLAGS="-Mpreprocess $FCFLAGS"
|
||||
fi
|
||||
|
||||
## Enable GPU offloading
|
||||
|
||||
# GPU offloading
|
||||
AC_ARG_ENABLE(gpu, [AS_HELP_STRING([--enable-gpu],[openmp|openacc : Use GPU-offloaded functions])], enable_gpu=$enableval, enable_gpu=no)
|
||||
AS_IF([test "$enable_gpu" = "yes"], [enable_gpu="openmp"])
|
||||
|
||||
# OpenMP offloading
|
||||
HAVE_OPENMP_OFFLOAD="no"
|
||||
AS_IF([test "$enable_gpu" = "openmp"], [
|
||||
AC_DEFINE([HAVE_OPENMP_OFFLOAD], [1], [If defined, activate OpenMP-offloaded routines])
|
||||
HAVE_OPENMP_OFFLOAD="yes"
|
||||
case $CC in
|
||||
|
||||
*gcc*)
|
||||
CFLAGS="$CFLAGS -fopenmp"
|
||||
;;
|
||||
*nvc*)
|
||||
CFLAGS="$CFLAGS -mp=gpu"
|
||||
;;
|
||||
esac
|
||||
|
||||
case $FC in
|
||||
|
||||
*gfortran*)
|
||||
FCFLAGS="$FCFLAGS -fopenmp"
|
||||
;;
|
||||
*nvfortran*)
|
||||
FCFLAGS="$FCFLAGS -mp=gpu"
|
||||
;;
|
||||
esac]
|
||||
)
|
||||
|
||||
# OpenMP offloading
|
||||
HAVE_OPENACC_OFFLOAD="no"
|
||||
AS_IF([test "$enable_gpu" = "openacc"], [
|
||||
AC_DEFINE([HAVE_OPENACC_OFFLOAD], [1], [If defined, activate OpenACC-offloaded routines])
|
||||
HAVE_OPENACC_OFFLOAD="yes"
|
||||
case $CC in
|
||||
|
||||
*gcc*)
|
||||
CFLAGS="$CFLAGS -fopenacc"
|
||||
;;
|
||||
*nvc*)
|
||||
CFLAGS="$CFLAGS -acc=gpu"
|
||||
;;
|
||||
esac
|
||||
|
||||
case $FC in
|
||||
|
||||
*gfortran*)
|
||||
FCFLAGS="$FCFLAGS -fopenacc"
|
||||
;;
|
||||
*nvfortran*)
|
||||
FCFLAGS="$FCFLAGS -acc=gpu"
|
||||
;;
|
||||
esac
|
||||
|
||||
])
|
||||
|
||||
# cuBLAS offloading
|
||||
AC_ARG_WITH(cublas, [AS_HELP_STRING([--with-cublas],[Use cuBLAS-offloaded functions])], HAVE_CUBLAS_OFFLOAD=$withval, HAVE_CUBLAS_OFFLOAD=no)
|
||||
AS_IF([test "$HAVE_CUBLAS_OFFLOAD" = "yes"], [
|
||||
AC_DEFINE([HAVE_CUBLAS_OFFLOAD], [1], [If defined, activate cuBLAS-offloaded routines])
|
||||
HAVE_OPENACC_OFFLOAD="yes"
|
||||
case $CC in
|
||||
|
||||
*gcc*)
|
||||
CFLAGS="$CFLAGS -fopenmp"
|
||||
LDFLAGS="-lcublas"
|
||||
;;
|
||||
*nvc*)
|
||||
CFLAGS="$CFLAGS -mp=gpu -cudalib=cublas"
|
||||
;;
|
||||
esac
|
||||
|
||||
case $FC in
|
||||
|
||||
*gfortran*)
|
||||
FCFLAGS="$FCFLAGS -fopenmp"
|
||||
;;
|
||||
*nvfortran*)
|
||||
FCFLAGS="$FCFLAGS -mp=gpu -cudalib=cublas"
|
||||
;;
|
||||
esac
|
||||
])
|
||||
|
||||
AC_ARG_ENABLE(malloc-trace, [AS_HELP_STRING([--enable-malloc-trace],[use debug malloc/free])], ok=$enableval, ok=no)
|
||||
if test "$ok" = "yes"; then
|
||||
AC_DEFINE(MALLOC_TRACE,"malloc_trace.dat",[Define to use debugging malloc/free])
|
||||
ARGS="${ARGS} malloc-trace"
|
||||
fi
|
||||
|
||||
AC_ARG_ENABLE(prof, [AS_HELP_STRING([--enable-prof],[compile for profiling])], ok=$enableval, ok=no)
|
||||
if test "$ok" = "yes"; then
|
||||
CFLAGS="${CFLAGS} -pg"
|
||||
AC_DEFINE(ENABLE_PROF,1,[Define when using the profiler tool])
|
||||
ARGS="${ARGS} prof"
|
||||
fi
|
||||
|
||||
AC_ARG_WITH(efence, [AS_HELP_STRING([--with-efence],[use ElectricFence library])], ok=$withval, ok=no)
|
||||
if test "$ok" = "yes"; then
|
||||
AC_CHECK_LIB([efence], [malloc])
|
||||
ARGS="${ARGS} efence"
|
||||
fi
|
||||
|
||||
|
||||
|
||||
##
|
||||
|
||||
AC_ARG_ENABLE(debug, [AS_HELP_STRING([--enable-debug],[compile for debugging])], ok=$enableval, ok=no)
|
||||
if test "$ok" = "yes"; then
|
||||
if test "$GCC" = "yes"; then
|
||||
@ -247,26 +376,6 @@ if test "$ok" = "yes"; then
|
||||
ARGS="${ARGS} debug"
|
||||
fi
|
||||
|
||||
AC_ARG_ENABLE(malloc-trace, [AS_HELP_STRING([--enable-malloc-trace],[use debug malloc/free])], ok=$enableval, ok=no)
|
||||
if test "$ok" = "yes"; then
|
||||
AC_DEFINE(MALLOC_TRACE,"malloc_trace.dat",[Define to use debugging malloc/free])
|
||||
ARGS="${ARGS} malloc-trace"
|
||||
fi
|
||||
|
||||
AC_ARG_ENABLE(prof, [AS_HELP_STRING([--enable-prof],[compile for profiling])], ok=$enableval, ok=no)
|
||||
if test "$ok" = "yes"; then
|
||||
CFLAGS="${CFLAGS} -pg"
|
||||
AC_DEFINE(ENABLE_PROF,1,[Define when using the profiler tool])
|
||||
ARGS="${ARGS} prof"
|
||||
fi
|
||||
|
||||
AC_ARG_WITH(efence, [AS_HELP_STRING([--with-efence],[use ElectricFence library])], ok=$withval, ok=no)
|
||||
if test "$ok" = "yes"; then
|
||||
AC_CHECK_LIB(efence, malloc)
|
||||
ARGS="${ARGS} efence"
|
||||
fi
|
||||
|
||||
|
||||
# Checks for header files.
|
||||
|
||||
## qmckl
|
||||
@ -321,21 +430,6 @@ if test "x${QMCKL_DEVEL}" != "x"; then
|
||||
|
||||
fi
|
||||
|
||||
# Enable Verificarlo tests
|
||||
AC_ARG_ENABLE([vfc_ci],
|
||||
[ --enable-vfc_ci Build the library with vfc_ci support],
|
||||
[case "${enableval}" in
|
||||
yes) vfc_ci=true && FCFLAGS="-D VFC_CI $FCFLAGS" && CFLAGS="-D VFC_CI $CFLAGS";;
|
||||
no) vfc_ci=false ;;
|
||||
*) AC_MSG_ERROR([bad value ${enableval} for --enable_vfc_ci]) ;;
|
||||
esac],[vfc_ci=false])
|
||||
AM_CONDITIONAL([VFC_CI], [test x$vfc_ci = xtrue])
|
||||
|
||||
if test "$FC" = "verificarlo-f"; then
|
||||
AC_MSG_NOTICE(verificarlo-f detected)
|
||||
# Arguments order is important here
|
||||
FCFLAGS="-Mpreprocess $FCFLAGS"
|
||||
fi
|
||||
|
||||
#PKG-CONFIG
|
||||
#mkl-dynamic-lp64-seq
|
||||
@ -371,6 +465,9 @@ LDFLAGS:........: ${LDFLAGS}
|
||||
LIBS............: ${LIBS}
|
||||
USE CHAMELEON...: ${with_chameleon}
|
||||
HPC version.....: ${HAVE_HPC}
|
||||
OpenMP offload..: ${HAVE_OPENMP_OFFLOAD}
|
||||
OpenACC offload.: ${HAVE_OPENACC_OFFLOAD}
|
||||
cuBLAS offload..: ${HAVE_CUBLAS_OFFLOAD}
|
||||
|
||||
Package features:
|
||||
${ARGS}
|
||||
|
199
org/examples.org
199
org/examples.org
@ -1,199 +0,0 @@
|
||||
#+TITLE: Code examples
|
||||
#+SETUPFILE: ../tools/theme.setup
|
||||
#+INCLUDE: ../tools/lib.org
|
||||
|
||||
In this section, we present examples of usage of QMCkl.
|
||||
For simplicity, we assume that the wave function parameters are stores
|
||||
in a [[https://github.com/TREX-CoE/trexio][TREXIO]] file.
|
||||
|
||||
* Checking errors
|
||||
|
||||
All QMCkl functions return an error code. A convenient way to handle
|
||||
errors is to write an error-checking function that displays the
|
||||
error in text format and exits the program.
|
||||
|
||||
#+NAME: qmckl_check_error
|
||||
#+begin_src f90
|
||||
subroutine qmckl_check_error(rc, message)
|
||||
use qmckl
|
||||
implicit none
|
||||
integer(qmckl_exit_code), intent(in) :: rc
|
||||
character(len=*) , intent(in) :: message
|
||||
character(len=128) :: str_buffer
|
||||
if (rc /= QMCKL_SUCCESS) then
|
||||
print *, message
|
||||
call qmckl_string_of_error(rc, str_buffer)
|
||||
print *, str_buffer
|
||||
call exit(rc)
|
||||
end if
|
||||
end subroutine qmckl_check_error
|
||||
#+end_src
|
||||
|
||||
* Computing an atomic orbital on a grid
|
||||
:PROPERTIES:
|
||||
:header-args: :tangle ao_grid.f90
|
||||
:END:
|
||||
|
||||
The following program, in Fortran, computes the values of an atomic
|
||||
orbital on a regular 3-dimensional grid. The 100^3 grid points are
|
||||
automatically defined, such that the molecule fits in a box with 5
|
||||
atomic units in the borders.
|
||||
|
||||
This program uses the ~qmckl_check_error~ function defined above.
|
||||
|
||||
To use this program, run
|
||||
|
||||
#+begin_src bash :tangle no
|
||||
$ ao_grid <trexio_file> <AO_id> <point_num>
|
||||
#+end_src
|
||||
|
||||
|
||||
#+begin_src f90 :noweb yes
|
||||
<<qmckl_check_error>>
|
||||
|
||||
program ao_grid
|
||||
use qmckl
|
||||
implicit none
|
||||
|
||||
integer(qmckl_context) :: qmckl_ctx ! QMCkl context
|
||||
integer(qmckl_exit_code) :: rc ! Exit code of QMCkl functions
|
||||
|
||||
character(len=128) :: trexio_filename
|
||||
character(len=128) :: str_buffer
|
||||
integer :: ao_id
|
||||
integer :: point_num_x
|
||||
|
||||
integer(c_int64_t) :: nucl_num
|
||||
double precision, allocatable :: nucl_coord(:,:)
|
||||
|
||||
integer(c_int64_t) :: point_num
|
||||
integer(c_int64_t) :: ao_num
|
||||
integer(c_int64_t) :: ipoint, i, j, k
|
||||
double precision :: x, y, z, dr(3)
|
||||
double precision :: rmin(3), rmax(3)
|
||||
double precision, allocatable :: points(:,:)
|
||||
double precision, allocatable :: ao_vgl(:,:,:)
|
||||
#+end_src
|
||||
|
||||
Start by fetching the command-line arguments:
|
||||
|
||||
#+begin_src f90
|
||||
if (iargc() /= 3) then
|
||||
print *, 'Syntax: ao_grid <trexio_file> <AO_id> <point_num>'
|
||||
call exit(-1)
|
||||
end if
|
||||
call getarg(1, trexio_filename)
|
||||
call getarg(2, str_buffer)
|
||||
read(str_buffer, *) ao_id
|
||||
call getarg(3, str_buffer)
|
||||
read(str_buffer, *) point_num_x
|
||||
|
||||
if (point_num_x < 0 .or. point_num_x > 300) then
|
||||
print *, 'Error: 0 < point_num < 300'
|
||||
call exit(-1)
|
||||
end if
|
||||
#+end_src
|
||||
|
||||
Create the QMCkl context and initialize it with the wave function
|
||||
present in the TREXIO file:
|
||||
|
||||
#+begin_src f90
|
||||
qmckl_ctx = qmckl_context_create()
|
||||
rc = qmckl_trexio_read(qmckl_ctx, trexio_filename, 1_8*len(trim(trexio_filename)))
|
||||
call qmckl_check_error(rc, 'Read TREXIO')
|
||||
#+end_src
|
||||
|
||||
We need to check that ~ao_id~ is in the range, so we get the total
|
||||
number of AOs from QMCkl:
|
||||
|
||||
#+begin_src f90
|
||||
rc = qmckl_get_ao_basis_ao_num(qmckl_ctx, ao_num)
|
||||
call qmckl_check_error(rc, 'Getting ao_num')
|
||||
|
||||
if (ao_id < 0 .or. ao_id > ao_num) then
|
||||
print *, 'Error: 0 < ao_id < ', ao_num
|
||||
call exit(-1)
|
||||
end if
|
||||
#+end_src
|
||||
|
||||
Now we will compute the limits of the box in which the molecule fits.
|
||||
For that, we first need to ask QMCkl the coordinates of nuclei.
|
||||
|
||||
#+begin_src f90
|
||||
rc = qmckl_get_nucleus_num(qmckl_ctx, nucl_num)
|
||||
call qmckl_check_error(rc, 'Get nucleus num')
|
||||
|
||||
allocate( nucl_coord(3, nucl_num) )
|
||||
rc = qmckl_get_nucleus_coord(qmckl_ctx, 'N', nucl_coord, 3_8*nucl_num)
|
||||
call qmckl_check_error(rc, 'Get nucleus coord')
|
||||
#+end_src
|
||||
|
||||
We now compute the coordinates of opposite points of the box, and
|
||||
the distance between points along the 3 directions:
|
||||
|
||||
#+begin_src f90
|
||||
rmin(1) = minval( nucl_coord(1,:) ) - 5.d0
|
||||
rmin(2) = minval( nucl_coord(2,:) ) - 5.d0
|
||||
rmin(3) = minval( nucl_coord(3,:) ) - 5.d0
|
||||
|
||||
rmax(1) = maxval( nucl_coord(1,:) ) + 5.d0
|
||||
rmax(2) = maxval( nucl_coord(2,:) ) + 5.d0
|
||||
rmax(3) = maxval( nucl_coord(3,:) ) + 5.d0
|
||||
|
||||
dr(1:3) = (rmax(1:3) - rmin(1:3)) / dble(point_num_x-1)
|
||||
#+end_src
|
||||
|
||||
We now produce the list of point coordinates where the AO will be
|
||||
evaluated:
|
||||
|
||||
#+begin_src f90
|
||||
point_num = point_num_x**3
|
||||
allocate( points(point_num, 3) )
|
||||
ipoint=0
|
||||
z = rmin(3)
|
||||
do k=1,point_num_x
|
||||
y = rmin(2)
|
||||
do j=1,point_num_x
|
||||
x = rmin(1)
|
||||
do i=1,point_num_x
|
||||
ipoint = ipoint+1
|
||||
points(ipoint,1) = x
|
||||
points(ipoint,2) = y
|
||||
points(ipoint,3) = z
|
||||
x = x + dr(1)
|
||||
end do
|
||||
y = y + dr(2)
|
||||
end do
|
||||
z = z + dr(3)
|
||||
end do
|
||||
#+end_src
|
||||
|
||||
We give the points to QMCkl:
|
||||
|
||||
#+begin_src f90
|
||||
rc = qmckl_set_point(qmckl_ctx, 'T', points, point_num)
|
||||
call qmckl_check_error(rc, 'Setting points')
|
||||
#+end_src
|
||||
|
||||
We allocate the space required to retrieve the values, gradients and
|
||||
Laplacian of all AOs, and ask to retrieve the values of the
|
||||
AOs computed at the point positions.
|
||||
|
||||
#+begin_src f90
|
||||
allocate( ao_vgl(ao_num, 5, point_num) )
|
||||
rc = qmckl_get_ao_basis_ao_vgl(qmckl_ctx, ao_vgl, ao_num*5_8*point_num)
|
||||
call qmckl_check_error(rc, 'Setting points')
|
||||
#+end_src
|
||||
|
||||
We finally print the value of the AO:
|
||||
|
||||
#+begin_src f90
|
||||
do ipoint=1, point_num
|
||||
print '(3(F16.10,X),E20.10)', points(ipoint, 1:3), ao_vgl(ao_id,1,ipoint)
|
||||
end do
|
||||
#+end_src
|
||||
|
||||
#+begin_src f90
|
||||
deallocate( nucl_coord, points, ao_vgl )
|
||||
end program ao_grid
|
||||
#+end_src
|
1536
org/qmckl_ao.org
1536
org/qmckl_ao.org
File diff suppressed because it is too large
Load Diff
@ -84,8 +84,8 @@ are not intended to be passed to external codes.
|
||||
|
||||
#+begin_src c :comments org :tangle (eval h_private_type) :exports none
|
||||
typedef struct qmckl_vector {
|
||||
int64_t size;
|
||||
double* restrict data;
|
||||
int64_t size;
|
||||
} qmckl_vector;
|
||||
#+end_src
|
||||
|
||||
@ -132,6 +132,13 @@ qmckl_exit_code
|
||||
qmckl_vector_free( qmckl_context context,
|
||||
qmckl_vector* vector)
|
||||
{
|
||||
if (vector == NULL) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_INVALID_ARG_2,
|
||||
"qmckl_vector_free",
|
||||
"Null pointer");
|
||||
}
|
||||
|
||||
/* Always true */
|
||||
assert (vector->data != NULL);
|
||||
|
||||
@ -160,8 +167,8 @@ qmckl_vector_free( qmckl_context context,
|
||||
|
||||
#+begin_src c :comments org :tangle (eval h_private_type) :exports none
|
||||
typedef struct qmckl_matrix {
|
||||
int64_t size[2];
|
||||
double* restrict data;
|
||||
int64_t size[2];
|
||||
} qmckl_matrix;
|
||||
#+end_src
|
||||
|
||||
@ -213,6 +220,13 @@ qmckl_exit_code
|
||||
qmckl_matrix_free( qmckl_context context,
|
||||
qmckl_matrix* matrix)
|
||||
{
|
||||
if (matrix == NULL) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_INVALID_ARG_2,
|
||||
"qmckl_matrix_free",
|
||||
"Null pointer");
|
||||
}
|
||||
|
||||
/* Always true */
|
||||
assert (matrix->data != NULL);
|
||||
|
||||
@ -245,9 +259,9 @@ qmckl_matrix_free( qmckl_context context,
|
||||
#define QMCKL_TENSOR_ORDER_MAX 16
|
||||
|
||||
typedef struct qmckl_tensor {
|
||||
double* restrict data;
|
||||
int64_t order;
|
||||
int64_t size[QMCKL_TENSOR_ORDER_MAX];
|
||||
double* restrict data;
|
||||
} qmckl_tensor;
|
||||
#+end_src
|
||||
|
||||
@ -307,6 +321,13 @@ qmckl_exit_code
|
||||
qmckl_tensor_free( qmckl_context context,
|
||||
qmckl_tensor* tensor)
|
||||
{
|
||||
if (tensor == NULL) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_INVALID_ARG_2,
|
||||
"qmckl_tensor_free",
|
||||
"Null pointer");
|
||||
}
|
||||
|
||||
/* Always true */
|
||||
assert (tensor->data != NULL);
|
||||
|
||||
|
@ -169,7 +169,7 @@ qmckl_context qmckl_context_check(const qmckl_context context) {
|
||||
if (context == QMCKL_NULL_CONTEXT)
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
|
||||
const qmckl_context_struct* const ctx = (const qmckl_context_struct*) context;
|
||||
const qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
|
||||
/* Try to access memory */
|
||||
if (ctx->tag != VALID_TAG) {
|
||||
@ -267,7 +267,7 @@ qmckl_context qmckl_context_create() {
|
||||
{
|
||||
ctx->tag = VALID_TAG;
|
||||
|
||||
const qmckl_context context = (const qmckl_context) ctx;
|
||||
const qmckl_context context = (qmckl_context) ctx;
|
||||
assert ( qmckl_context_check(context) != QMCKL_NULL_CONTEXT );
|
||||
|
||||
qmckl_exit_code rc;
|
||||
|
@ -182,7 +182,7 @@ qmckl_exit_code qmckl_init_determinant(qmckl_context context) {
|
||||
return false;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
ctx->det.uninitialized = (1 << 6) - 1;
|
||||
@ -216,7 +216,7 @@ bool qmckl_determinant_provided(const qmckl_context context) {
|
||||
return false;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
return ctx->det.provided;
|
||||
@ -238,7 +238,7 @@ char qmckl_get_determinant_type (const qmckl_context context) {
|
||||
return (char) 0;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1;
|
||||
@ -256,7 +256,7 @@ int64_t qmckl_get_determinant_walk_num (const qmckl_context context) {
|
||||
return (int64_t) 0;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 1;
|
||||
@ -274,7 +274,7 @@ int64_t qmckl_get_determinant_det_num_alpha (const qmckl_context context) {
|
||||
return (int64_t) 0;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 2;
|
||||
@ -292,7 +292,7 @@ int64_t qmckl_get_determinant_det_num_beta (const qmckl_context context) {
|
||||
return (int64_t) 0;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 3;
|
||||
@ -310,7 +310,7 @@ int64_t* qmckl_get_determinant_mo_index_alpha (const qmckl_context context) {
|
||||
return NULL;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 4;
|
||||
@ -328,7 +328,7 @@ int64_t* qmckl_get_determinant_mo_index_beta (const qmckl_context context) {
|
||||
return NULL;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 5;
|
||||
@ -363,7 +363,7 @@ if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
#+end_src
|
||||
|
||||
#+NAME:post2
|
||||
@ -525,7 +525,7 @@ qmckl_exit_code qmckl_finalize_determinant(qmckl_context context) {
|
||||
NULL);
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
qmckl_exit_code rc;
|
||||
@ -596,7 +596,7 @@ qmckl_exit_code qmckl_get_det_vgl_alpha(qmckl_context context, double * const de
|
||||
rc = qmckl_provide_det_vgl_alpha(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = 5 * ctx->det.det_num_alpha * ctx->det.walk_num *
|
||||
@ -623,7 +623,7 @@ qmckl_exit_code qmckl_get_det_vgl_beta(qmckl_context context, double * const det
|
||||
rc = qmckl_provide_det_vgl_beta(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = 5 * ctx->det.det_num_beta * ctx->det.walk_num *
|
||||
@ -649,7 +649,7 @@ qmckl_exit_code qmckl_provide_det_vgl_alpha(qmckl_context context) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if(!(ctx->nucleus.provided)) {
|
||||
@ -748,7 +748,7 @@ qmckl_exit_code qmckl_provide_det_vgl_beta(qmckl_context context) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if(!(ctx->nucleus.provided)) {
|
||||
@ -1134,36 +1134,28 @@ end function qmckl_compute_det_vgl_beta_f
|
||||
|
||||
#+begin_src c :tangle (eval c_test) :exports none
|
||||
|
||||
#define walk_num chbrclf_walk_num
|
||||
#define elec_num chbrclf_elec_num
|
||||
#define shell_num chbrclf_shell_num
|
||||
#define ao_num chbrclf_ao_num
|
||||
|
||||
int64_t elec_up_num = chbrclf_elec_up_num;
|
||||
int64_t elec_dn_num = chbrclf_elec_dn_num;
|
||||
double* elec_coord = &(chbrclf_elec_coord[0][0][0]);
|
||||
const int64_t nucl_num = chbrclf_nucl_num;
|
||||
const double* nucl_charge = chbrclf_charge;
|
||||
const double* nucl_coord = &(chbrclf_nucl_coord[0][0]);
|
||||
|
||||
rc = qmckl_set_electron_num (context, elec_up_num, elec_dn_num);
|
||||
rc = qmckl_set_electron_num (context, chbrclf_elec_up_num, chbrclf_elec_dn_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_electron_walk_num (context, walk_num);
|
||||
rc = qmckl_set_electron_walk_num (context, chbrclf_walk_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
assert(qmckl_electron_provided(context));
|
||||
|
||||
rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*elec_num*3);
|
||||
rc = qmckl_set_electron_coord (context, 'N', elec_coord, chbrclf_walk_num*chbrclf_elec_num*3);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_num (context, nucl_num);
|
||||
rc = qmckl_set_nucleus_num (context, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0]), nucl_num*3);
|
||||
rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0]), chbrclf_nucl_num*3);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_charge(context, nucl_charge, nucl_num);
|
||||
rc = qmckl_set_nucleus_charge(context, nucl_charge, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
assert(qmckl_nucleus_provided(context));
|
||||
@ -1195,27 +1187,27 @@ rc = qmckl_set_ao_basis_prim_num (context, chbrclf_prim_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_ao_basis_provided(context));
|
||||
|
||||
rc = qmckl_set_ao_basis_nucleus_index (context, nucleus_index, nucl_num);
|
||||
rc = qmckl_set_ao_basis_nucleus_index (context, nucleus_index, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_ao_basis_provided(context));
|
||||
|
||||
rc = qmckl_set_ao_basis_nucleus_shell_num (context, nucleus_shell_num, nucl_num);
|
||||
rc = qmckl_set_ao_basis_nucleus_shell_num (context, nucleus_shell_num, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_ao_basis_provided(context));
|
||||
|
||||
rc = qmckl_set_ao_basis_shell_ang_mom (context, shell_ang_mom, shell_num);
|
||||
rc = qmckl_set_ao_basis_shell_ang_mom (context, shell_ang_mom, chbrclf_shell_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_ao_basis_provided(context));
|
||||
|
||||
rc = qmckl_set_ao_basis_shell_factor (context, shell_factor, shell_num);
|
||||
rc = qmckl_set_ao_basis_shell_factor (context, shell_factor, chbrclf_shell_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_ao_basis_provided(context));
|
||||
|
||||
rc = qmckl_set_ao_basis_shell_prim_num (context, shell_prim_num, shell_num);
|
||||
rc = qmckl_set_ao_basis_shell_prim_num (context, shell_prim_num, chbrclf_shell_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_ao_basis_provided(context));
|
||||
|
||||
rc = qmckl_set_ao_basis_shell_prim_index (context, shell_prim_index, shell_num);
|
||||
rc = qmckl_set_ao_basis_shell_prim_index (context, shell_prim_index, chbrclf_shell_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_ao_basis_provided(context));
|
||||
|
||||
@ -1239,14 +1231,13 @@ assert(rc == QMCKL_SUCCESS);
|
||||
assert(qmckl_ao_basis_provided(context));
|
||||
|
||||
|
||||
double ao_vgl[walk_num*elec_num][5][chbrclf_ao_num];
|
||||
double ao_vgl[chbrclf_walk_num*chbrclf_elec_num][5][chbrclf_ao_num];
|
||||
|
||||
rc = qmckl_get_ao_basis_ao_vgl(context, &(ao_vgl[0][0][0]), (int64_t) 5*walk_num*elec_num*chbrclf_ao_num);
|
||||
rc = qmckl_get_ao_basis_ao_vgl(context, &(ao_vgl[0][0][0]), (int64_t) 5*chbrclf_walk_num*chbrclf_elec_num*chbrclf_ao_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
/* Set up MO data */
|
||||
const int64_t mo_num = chbrclf_mo_num;
|
||||
rc = qmckl_set_mo_basis_mo_num(context, mo_num);
|
||||
rc = qmckl_set_mo_basis_mo_num(context, chbrclf_mo_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
const double * mo_coefficient = &(chbrclf_mo_coef[0]);
|
||||
@ -1256,31 +1247,31 @@ assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
assert(qmckl_mo_basis_provided(context));
|
||||
|
||||
double mo_vgl[walk_num*elec_num][5][chbrclf_mo_num];
|
||||
rc = qmckl_get_mo_basis_mo_vgl(context, &(mo_vgl[0][0][0]), 5*walk_num*elec_num*chbrclf_mo_num);
|
||||
double mo_vgl[chbrclf_walk_num*chbrclf_elec_num][5][chbrclf_mo_num];
|
||||
rc = qmckl_get_mo_basis_mo_vgl(context, &(mo_vgl[0][0][0]), 5*chbrclf_walk_num*chbrclf_elec_num*chbrclf_mo_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
/* Set up determinant data */
|
||||
|
||||
const int64_t det_num_alpha = 1;
|
||||
const int64_t det_num_beta = 1;
|
||||
int64_t mo_index_alpha[det_num_alpha][walk_num][elec_up_num];
|
||||
int64_t mo_index_beta[det_num_alpha][walk_num][elec_dn_num];
|
||||
#define det_num_alpha 1
|
||||
#define det_num_beta 1
|
||||
int64_t mo_index_alpha[det_num_alpha][chbrclf_walk_num][chbrclf_elec_up_num];
|
||||
int64_t mo_index_beta[det_num_alpha][chbrclf_walk_num][chbrclf_elec_dn_num];
|
||||
|
||||
int i, j, k;
|
||||
for(k = 0; k < det_num_alpha; ++k)
|
||||
for(i = 0; i < walk_num; ++i)
|
||||
for(j = 0; j < elec_up_num; ++j)
|
||||
for(i = 0; i < chbrclf_walk_num; ++i)
|
||||
for(j = 0; j < chbrclf_elec_up_num; ++j)
|
||||
mo_index_alpha[k][i][j] = j + 1;
|
||||
for(k = 0; k < det_num_beta; ++k)
|
||||
for(i = 0; i < walk_num; ++i)
|
||||
for(j = 0; j < elec_up_num; ++j)
|
||||
for(i = 0; i < chbrclf_walk_num; ++i)
|
||||
for(j = 0; j < chbrclf_elec_up_num; ++j)
|
||||
mo_index_beta[k][i][j] = j + 1;
|
||||
|
||||
rc = qmckl_set_determinant_type (context, typ);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_determinant_walk_num (context, walk_num);
|
||||
rc = qmckl_set_determinant_walk_num (context, chbrclf_walk_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_determinant_det_num_alpha (context, det_num_alpha);
|
||||
@ -1297,8 +1288,8 @@ assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
// Get slater-determinant
|
||||
|
||||
double det_vgl_alpha[det_num_alpha][walk_num][5][elec_up_num][elec_up_num];
|
||||
double det_vgl_beta[det_num_beta][walk_num][5][elec_dn_num][elec_dn_num];
|
||||
double det_vgl_alpha[det_num_alpha][chbrclf_walk_num][5][chbrclf_elec_up_num][chbrclf_elec_up_num];
|
||||
double det_vgl_beta[det_num_beta][chbrclf_walk_num][5][chbrclf_elec_dn_num][chbrclf_elec_dn_num];
|
||||
|
||||
rc = qmckl_get_det_vgl_alpha(context, &(det_vgl_alpha[0][0][0][0][0]));
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
@ -1347,7 +1338,7 @@ qmckl_exit_code qmckl_get_det_inv_matrix_alpha(qmckl_context context, double * c
|
||||
rc = qmckl_provide_det_inv_matrix_alpha(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->det.det_num_alpha * ctx->det.walk_num * ctx->electron.up_num * ctx->electron.up_num;
|
||||
@ -1376,7 +1367,7 @@ qmckl_exit_code qmckl_get_det_inv_matrix_beta(qmckl_context context, double * co
|
||||
rc = qmckl_provide_det_inv_matrix_beta(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->det.det_num_alpha * ctx->det.walk_num * ctx->electron.down_num * ctx->electron.down_num;
|
||||
@ -1405,7 +1396,7 @@ qmckl_exit_code qmckl_get_det_adj_matrix_alpha(qmckl_context context, double * c
|
||||
rc = qmckl_provide_det_inv_matrix_alpha(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->det.det_num_alpha * ctx->det.walk_num * ctx->electron.up_num * ctx->electron.up_num;
|
||||
@ -1434,7 +1425,7 @@ qmckl_exit_code qmckl_get_det_adj_matrix_beta(qmckl_context context, double * co
|
||||
rc = qmckl_provide_det_inv_matrix_beta(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->det.det_num_alpha * ctx->det.walk_num * ctx->electron.down_num * ctx->electron.down_num;
|
||||
@ -1463,7 +1454,7 @@ qmckl_exit_code qmckl_get_det_alpha(qmckl_context context, double * const det_va
|
||||
rc = qmckl_provide_det_inv_matrix_alpha(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->det.det_num_alpha * ctx->det.walk_num;
|
||||
@ -1492,7 +1483,7 @@ qmckl_exit_code qmckl_get_det_beta(qmckl_context context, double * const det_val
|
||||
rc = qmckl_provide_det_inv_matrix_beta(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->det.det_num_alpha * ctx->det.walk_num;
|
||||
@ -1517,7 +1508,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_alpha(qmckl_context context) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if(!(ctx->nucleus.provided)) {
|
||||
@ -1640,7 +1631,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_beta(qmckl_context context) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if(!(ctx->nucleus.provided)) {
|
||||
@ -2047,8 +2038,8 @@ end function qmckl_compute_det_inv_matrix_beta_f
|
||||
#+begin_src c :tangle (eval c_test) :exports none
|
||||
// Get adjoint of the slater-determinant
|
||||
|
||||
double det_inv_matrix_alpha[det_num_alpha][walk_num][elec_up_num][elec_up_num];
|
||||
double det_inv_matrix_beta[det_num_beta][walk_num][elec_dn_num][elec_dn_num];
|
||||
double det_inv_matrix_alpha[det_num_alpha][chbrclf_walk_num][chbrclf_elec_up_num][chbrclf_elec_up_num];
|
||||
double det_inv_matrix_beta[det_num_beta][chbrclf_walk_num][chbrclf_elec_dn_num][chbrclf_elec_dn_num];
|
||||
|
||||
rc = qmckl_get_det_inv_matrix_alpha(context, &(det_inv_matrix_alpha[0][0][0][0]));
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
@ -133,7 +133,7 @@ integer function qmckl_distance_sq_f(context, transa, transb, m, n, &
|
||||
|
||||
if (transb == 'N' .or. transb == 'n') then
|
||||
continue
|
||||
else if (transa == 'T' .or. transa == 't') then
|
||||
else if (transb == 'T' .or. transb == 't') then
|
||||
transab = transab + 2
|
||||
else
|
||||
transab = -100
|
||||
@ -533,7 +533,7 @@ integer function qmckl_distance_f(context, transa, transb, m, n, &
|
||||
|
||||
if (transb == 'N' .or. transb == 'n') then
|
||||
continue
|
||||
else if (transa == 'T' .or. transa == 't') then
|
||||
else if (transb == 'T' .or. transb == 't') then
|
||||
transab = transab + 2
|
||||
else
|
||||
transab = -100
|
||||
@ -1314,7 +1314,7 @@ integer function qmckl_distance_rescaled_deriv_e_f(context, transa, transb, m, n
|
||||
|
||||
if (transb == 'N' .or. transb == 'n') then
|
||||
continue
|
||||
else if (transa == 'T' .or. transa == 't') then
|
||||
else if (transb == 'T' .or. transb == 't') then
|
||||
transab = transab + 2
|
||||
else
|
||||
transab = -100
|
||||
|
@ -97,8 +97,8 @@ int main() {
|
||||
| ~ee_distance_rescaled_deriv_e_date~ | ~uint64_t~ | Last modification date of the electron-electron distance derivatives |
|
||||
| ~ee_pot~ | ~double[walk_num]~ | Electron-electron rescaled distances derivatives |
|
||||
| ~ee_pot_date~ | ~uint64_t~ | Last modification date of the electron-electron distance derivatives |
|
||||
| ~en_pot~ | double[walk_num] | Electron-nucleus potential energy |
|
||||
| ~en_pot_date~ | int64_t | Date when the electron-nucleus potential energy was computed |
|
||||
| ~en_pot~ | ~double[walk_num]~ | Electron-nucleus potential energy |
|
||||
| ~en_pot_date~ | ~int64_t~ | Date when the electron-nucleus potential energy was computed |
|
||||
| ~en_distance_rescaled~ | ~double[walk_num][nucl_num][num]~ | Electron-nucleus distances |
|
||||
| ~en_distance_rescaled_date~ | ~uint64_t~ | Last modification date of the electron-electron distances |
|
||||
| ~en_distance_rescaled_deriv_e~ | ~double[walk_num][4][nucl_num][num]~ | Electron-electron rescaled distances derivatives |
|
||||
@ -157,7 +157,7 @@ qmckl_exit_code qmckl_init_electron(qmckl_context context) {
|
||||
return false;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
ctx->electron.uninitialized = (1 << 2) - 1;
|
||||
@ -182,7 +182,7 @@ bool qmckl_electron_provided(const qmckl_context context) {
|
||||
return false;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
return ctx->electron.provided;
|
||||
@ -228,7 +228,7 @@ qmckl_get_electron_num (const qmckl_context context, int64_t* const num) {
|
||||
"num is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 0;
|
||||
@ -256,7 +256,7 @@ qmckl_get_electron_up_num (const qmckl_context context, int64_t* const up_num) {
|
||||
"up_num is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 0;
|
||||
@ -284,7 +284,7 @@ qmckl_get_electron_down_num (const qmckl_context context, int64_t* const down_nu
|
||||
"down_num is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 0;
|
||||
@ -323,7 +323,7 @@ qmckl_get_electron_walk_num (const qmckl_context context, int64_t* const walk_nu
|
||||
"walk_num is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 1;
|
||||
@ -360,7 +360,7 @@ qmckl_get_electron_rescale_factor_ee (const qmckl_context context, double* const
|
||||
"rescale_factor_kappa_ee is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
assert (ctx->electron.rescale_factor_kappa_ee > 0.0);
|
||||
@ -383,7 +383,7 @@ qmckl_get_electron_rescale_factor_en (const qmckl_context context, double* const
|
||||
"rescale_factor_kappa_en is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
assert (ctx->electron.rescale_factor_kappa_en > 0.0);
|
||||
@ -448,7 +448,7 @@ qmckl_get_electron_coord (const qmckl_context context,
|
||||
return QMCKL_INVALID_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!ctx->electron.provided) {
|
||||
@ -489,7 +489,14 @@ if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
|
||||
if (mask != 0 && !(ctx->electron.uninitialized & mask)) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_ALREADY_SET,
|
||||
"qmckl_set_electron_*",
|
||||
NULL);
|
||||
}
|
||||
#+end_src
|
||||
|
||||
#+NAME:post2
|
||||
@ -544,6 +551,8 @@ qmckl_exit_code
|
||||
qmckl_set_electron_num(qmckl_context context,
|
||||
const int64_t up_num,
|
||||
const int64_t down_num) {
|
||||
int32_t mask = 1 << 0;
|
||||
|
||||
<<pre2>>
|
||||
|
||||
if (up_num <= 0) {
|
||||
@ -560,8 +569,6 @@ qmckl_set_electron_num(qmckl_context context,
|
||||
"down_num < 0");
|
||||
}
|
||||
|
||||
int32_t mask = 1 << 0;
|
||||
|
||||
ctx->electron.up_num = up_num;
|
||||
ctx->electron.down_num = down_num;
|
||||
ctx->electron.num = up_num + down_num;
|
||||
@ -576,6 +583,8 @@ qmckl_set_electron_num(qmckl_context context,
|
||||
qmckl_exit_code
|
||||
qmckl_set_electron_walk_num(qmckl_context context, const int64_t walk_num) {
|
||||
|
||||
int32_t mask = 1 << 1;
|
||||
|
||||
<<pre2>>
|
||||
|
||||
if (walk_num <= 0) {
|
||||
@ -585,7 +594,6 @@ qmckl_set_electron_walk_num(qmckl_context context, const int64_t walk_num) {
|
||||
"walk_num <= 0");
|
||||
}
|
||||
|
||||
int32_t mask = 1 << 1;
|
||||
ctx->electron.walk_num = walk_num;
|
||||
|
||||
<<post2>>
|
||||
@ -598,6 +606,9 @@ qmckl_set_electron_walk_num(qmckl_context context, const int64_t walk_num) {
|
||||
qmckl_exit_code
|
||||
qmckl_set_electron_rescale_factor_ee(qmckl_context context,
|
||||
const double rescale_factor_kappa_ee) {
|
||||
|
||||
int32_t mask = 0; // can be changed
|
||||
|
||||
<<pre2>>
|
||||
|
||||
if (rescale_factor_kappa_ee <= 0.0) {
|
||||
@ -615,6 +626,9 @@ qmckl_set_electron_rescale_factor_ee(qmckl_context context,
|
||||
qmckl_exit_code
|
||||
qmckl_set_electron_rescale_factor_en(qmckl_context context,
|
||||
const double rescale_factor_kappa_en) {
|
||||
|
||||
int32_t mask = 0; // can be changed
|
||||
|
||||
<<pre2>>
|
||||
|
||||
if (rescale_factor_kappa_en <= 0.0) {
|
||||
@ -675,6 +689,8 @@ qmckl_set_electron_coord(qmckl_context context,
|
||||
const int64_t size_max)
|
||||
{
|
||||
|
||||
int32_t mask = 0; // coord can be changed
|
||||
|
||||
<<pre2>>
|
||||
|
||||
if (transp != 'N' && transp != 'T') {
|
||||
@ -718,7 +734,7 @@ qmckl_set_electron_coord(qmckl_context context,
|
||||
ctx->electron.coord_old = ctx->electron.coord_new ;
|
||||
|
||||
qmckl_exit_code rc;
|
||||
rc = qmckl_set_point(context, transp, coord, size_max/3);
|
||||
rc = qmckl_set_point(context, transp, size_max/3, coord, size_max);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
ctx->electron.coord_new = ctx->point.coord ;
|
||||
@ -897,7 +913,7 @@ qmckl_exit_code qmckl_get_electron_ee_distance(qmckl_context context, double* co
|
||||
rc = qmckl_provide_ee_distance(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->electron.num * ctx->electron.num * ctx->electron.walk_num;
|
||||
@ -921,7 +937,7 @@ qmckl_exit_code qmckl_provide_ee_distance(qmckl_context context)
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
|
||||
@ -1138,7 +1154,7 @@ qmckl_exit_code qmckl_get_electron_ee_distance_rescaled(qmckl_context context, d
|
||||
rc = qmckl_provide_ee_distance_rescaled(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->electron.num * ctx->electron.num * ctx->electron.walk_num;
|
||||
@ -1162,7 +1178,7 @@ qmckl_exit_code qmckl_provide_ee_distance_rescaled(qmckl_context context)
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
|
||||
@ -1218,7 +1234,7 @@ qmckl_exit_code qmckl_provide_ee_distance_rescaled(qmckl_context context)
|
||||
| ~elec_num~ | ~int64_t~ | in | Number of electrons |
|
||||
| ~rescale_factor_kappa_ee~ | ~double~ | in | Factor to rescale ee distances |
|
||||
| ~walk_num~ | ~int64_t~ | in | Number of walkers |
|
||||
| ~coord~ | ~double[walk_num][3][elec_num]~ | in | Electron coordinates |
|
||||
| ~coord~ | ~double[3][walk_num][elec_num]~ | in | Electron coordinates |
|
||||
| ~ee_distance~ | ~double[walk_num][elec_num][elec_num]~ | out | Electron-electron rescaled distances |
|
||||
|
||||
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
|
||||
@ -1231,7 +1247,7 @@ integer function qmckl_compute_ee_distance_rescaled_f(context, elec_num, rescale
|
||||
integer*8 , intent(in) :: elec_num
|
||||
double precision , intent(in) :: rescale_factor_kappa_ee
|
||||
integer*8 , intent(in) :: walk_num
|
||||
double precision , intent(in) :: coord(elec_num,3,walk_num)
|
||||
double precision , intent(in) :: coord(elec_num,walk_num,3)
|
||||
double precision , intent(out) :: ee_distance_rescaled(elec_num,elec_num,walk_num)
|
||||
|
||||
integer*8 :: k
|
||||
@ -1357,7 +1373,7 @@ assert(fabs(ee_distance_rescaled[elec_num*elec_num+1]-0.9985724058042633) < 1.e-
|
||||
|
||||
#+end_src
|
||||
|
||||
** Electron-electron rescaled distance gradients and laplacian with respect to electron coords
|
||||
** Electron-electron rescaled distance gradients and Laplacian with respect to electron coords
|
||||
|
||||
The rescaled distances which is given as $R = (1 - \exp{-\kappa r})/\kappa$
|
||||
needs to be perturbed with respect to the electorn coordinates.
|
||||
@ -1384,7 +1400,7 @@ qmckl_exit_code qmckl_get_electron_ee_distance_rescaled_deriv_e(qmckl_context co
|
||||
rc = qmckl_provide_ee_distance_rescaled_deriv_e(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = 4 * ctx->electron.num * ctx->electron.num * ctx->electron.walk_num;
|
||||
@ -1408,7 +1424,7 @@ qmckl_exit_code qmckl_provide_ee_distance_rescaled_deriv_e(qmckl_context context
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
|
||||
@ -1464,7 +1480,7 @@ qmckl_exit_code qmckl_provide_ee_distance_rescaled_deriv_e(qmckl_context context
|
||||
| ~elec_num~ | ~int64_t~ | in | Number of electrons |
|
||||
| ~rescale_factor_kappa_ee~ | ~double~ | in | Factor to rescale ee distances |
|
||||
| ~walk_num~ | ~int64_t~ | in | Number of walkers |
|
||||
| ~coord~ | ~double[walk_num][3][elec_num]~ | in | Electron coordinates |
|
||||
| ~coord~ | ~double[3][walk_num][elec_num]~ | in | Electron coordinates |
|
||||
| ~ee_distance_deriv_e~ | ~double[walk_num][4][elec_num][elec_num]~ | out | Electron-electron rescaled distance derivatives |
|
||||
|
||||
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
|
||||
@ -1477,7 +1493,7 @@ integer function qmckl_compute_ee_distance_rescaled_deriv_e_f(context, elec_num,
|
||||
integer*8 , intent(in) :: elec_num
|
||||
double precision , intent(in) :: rescale_factor_kappa_ee
|
||||
integer*8 , intent(in) :: walk_num
|
||||
double precision , intent(in) :: coord(elec_num,3,walk_num)
|
||||
double precision , intent(in) :: coord(elec_num,walk_num,3)
|
||||
double precision , intent(out) :: ee_distance_rescaled_deriv_e(4,elec_num,elec_num,walk_num)
|
||||
|
||||
integer*8 :: k
|
||||
@ -1501,8 +1517,8 @@ integer function qmckl_compute_ee_distance_rescaled_deriv_e_f(context, elec_num,
|
||||
|
||||
do k=1,walk_num
|
||||
info = qmckl_distance_rescaled_deriv_e(context, 'T', 'T', elec_num, elec_num, &
|
||||
coord(1,1,k), elec_num, &
|
||||
coord(1,1,k), elec_num, &
|
||||
coord(1,k,1), elec_num*walk_num, &
|
||||
coord(1,k,1), elec_num*walk_num, &
|
||||
ee_distance_rescaled_deriv_e(1,1,1,k), elec_num, rescale_factor_kappa_ee)
|
||||
if (info /= QMCKL_SUCCESS) then
|
||||
exit
|
||||
@ -1613,7 +1629,7 @@ qmckl_exit_code qmckl_get_electron_ee_potential(qmckl_context context, double* c
|
||||
rc = qmckl_provide_ee_potential(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->electron.walk_num * sizeof(double);
|
||||
@ -1637,7 +1653,7 @@ qmckl_exit_code qmckl_provide_ee_potential(qmckl_context context)
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!ctx->electron.provided) return QMCKL_NOT_PROVIDED;
|
||||
@ -1818,7 +1834,7 @@ qmckl_exit_code qmckl_get_electron_en_distance(qmckl_context context, double* di
|
||||
rc = qmckl_provide_en_distance(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->electron.num * ctx->nucleus.num * ctx->electron.walk_num;
|
||||
@ -1842,7 +1858,7 @@ qmckl_exit_code qmckl_provide_en_distance(qmckl_context context)
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!(ctx->nucleus.provided)) {
|
||||
@ -1905,7 +1921,7 @@ qmckl_exit_code qmckl_provide_en_distance(qmckl_context context)
|
||||
| ~elec_num~ | ~int64_t~ | in | Number of electrons |
|
||||
| ~nucl_num~ | ~int64_t~ | in | Number of nuclei |
|
||||
| ~walk_num~ | ~int64_t~ | in | Number of walkers |
|
||||
| ~elec_coord~ | ~double[walk_num][3][elec_num]~ | in | Electron coordinates |
|
||||
| ~elec_coord~ | ~double[3][walk_num][elec_num]~ | in | Electron coordinates |
|
||||
| ~nucl_coord~ | ~double[3][elec_num]~ | in | Nuclear coordinates |
|
||||
| ~en_distance~ | ~double[walk_num][nucl_num][elec_num]~ | out | Electron-nucleus distances |
|
||||
|
||||
@ -2097,7 +2113,7 @@ qmckl_exit_code qmckl_get_electron_en_distance_rescaled(qmckl_context context, d
|
||||
rc = qmckl_provide_en_distance_rescaled(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->electron.num * ctx->nucleus.num * ctx->electron.walk_num;
|
||||
@ -2121,7 +2137,7 @@ qmckl_exit_code qmckl_provide_en_distance_rescaled(qmckl_context context)
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!(ctx->nucleus.provided)) {
|
||||
@ -2183,7 +2199,7 @@ qmckl_exit_code qmckl_provide_en_distance_rescaled(qmckl_context context)
|
||||
| ~nucl_num~ | ~int64_t~ | in | Number of nuclei |
|
||||
| ~rescale_factor_kappa_en~ | ~double~ | in | The factor for rescaled distances |
|
||||
| ~walk_num~ | ~int64_t~ | in | Number of walkers |
|
||||
| ~elec_coord~ | ~double[walk_num][3][elec_num]~ | in | Electron coordinates |
|
||||
| ~elec_coord~ | ~double[3][walk_num][elec_num]~ | in | Electron coordinates |
|
||||
| ~nucl_coord~ | ~double[3][elec_num]~ | in | Nuclear coordinates |
|
||||
| ~en_distance_rescaled~ | ~double[walk_num][nucl_num][elec_num]~ | out | Electron-nucleus distances |
|
||||
|
||||
@ -2318,16 +2334,6 @@ print ( "[1][0][1] : ", (1.0 - np.exp(-kappa * np.linalg.norm(elec_2_w2-nucl_1))
|
||||
#+begin_src c :tangle (eval c_test)
|
||||
|
||||
assert(qmckl_electron_provided(context));
|
||||
|
||||
rc = qmckl_set_nucleus_num (context, nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_charge (context, charge, nucl_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_coord (context, 'T', nucl_coord, 3*nucl_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
assert(qmckl_nucleus_provided(context));
|
||||
|
||||
double en_distance_rescaled[walk_num][nucl_num][elec_num];
|
||||
@ -2385,7 +2391,7 @@ qmckl_exit_code qmckl_get_electron_en_distance_rescaled_deriv_e(qmckl_context co
|
||||
rc = qmckl_provide_en_distance_rescaled_deriv_e(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = 4 * ctx->electron.num * ctx->nucleus.num * ctx->electron.walk_num;
|
||||
@ -2409,7 +2415,7 @@ qmckl_exit_code qmckl_provide_en_distance_rescaled_deriv_e(qmckl_context context
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!(ctx->nucleus.provided)) {
|
||||
@ -2471,9 +2477,9 @@ qmckl_exit_code qmckl_provide_en_distance_rescaled_deriv_e(qmckl_context context
|
||||
| ~nucl_num~ | ~int64_t~ | in | Number of nuclei |
|
||||
| ~rescale_factor_kappa_en~ | ~double~ | in | The factor for rescaled distances |
|
||||
| ~walk_num~ | ~int64_t~ | in | Number of walkers |
|
||||
| ~elec_coord~ | ~double[walk_num][3][elec_num]~ | in | Electron coordinates |
|
||||
| ~elec_coord~ | ~double[3][walk_num][elec_num]~ | in | Electron coordinates |
|
||||
| ~nucl_coord~ | ~double[3][elec_num]~ | in | Nuclear coordinates |
|
||||
| ~en_distance_rescaled_deriv_e~ | ~double[walk_num][4][nucl_num][elec_num]~ | out | Electron-nucleus distance derivatives |
|
||||
| ~en_distance_rescaled_deriv_e~ | ~double[walk_num][nucl_num][elec_num][4]~ | out | Electron-nucleus distance derivatives |
|
||||
|
||||
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
|
||||
integer function qmckl_compute_en_distance_rescaled_deriv_e_f(context, elec_num, nucl_num, &
|
||||
@ -2586,18 +2592,9 @@ import numpy as np
|
||||
|
||||
assert(qmckl_electron_provided(context));
|
||||
|
||||
rc = qmckl_set_nucleus_num (context, nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_rescale_factor (context, nucl_rescale_factor_kappa);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_charge (context, charge, nucl_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_coord (context, 'T', nucl_coord, 3*nucl_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
assert(qmckl_nucleus_provided(context));
|
||||
|
||||
double en_distance_rescaled_deriv_e[walk_num][4][nucl_num][elec_num];
|
||||
@ -2656,7 +2653,7 @@ qmckl_exit_code qmckl_get_electron_en_potential(qmckl_context context, double* c
|
||||
rc = qmckl_provide_en_potential(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->electron.walk_num * sizeof(double);
|
||||
@ -2680,7 +2677,7 @@ qmckl_exit_code qmckl_provide_en_potential(qmckl_context context)
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!ctx->electron.provided) return QMCKL_NOT_PROVIDED;
|
||||
|
@ -105,7 +105,8 @@ typedef int32_t qmckl_exit_code;
|
||||
| ~QMCKL_DEALLOCATION_FAILED~ | 105 | 'De-allocation failed' |
|
||||
| ~QMCKL_NOT_PROVIDED~ | 106 | 'Not provided' |
|
||||
| ~QMCKL_OUT_OF_BOUNDS~ | 107 | 'Index out of bounds' |
|
||||
| ~QMCKL_INVALID_EXIT_CODE~ | 108 | 'Invalid exit code' |
|
||||
| ~QMCKL_ALREADY_SET~ | 108 | 'Already set' |
|
||||
| ~QMCKL_INVALID_EXIT_CODE~ | 109 | 'Invalid exit code' |
|
||||
|
||||
# We need to force Emacs not to indent the Python code:
|
||||
# -*- org-src-preserve-indentation: t
|
||||
@ -164,7 +165,8 @@ return '\n'.join(result)
|
||||
#define QMCKL_DEALLOCATION_FAILED ((qmckl_exit_code) 105)
|
||||
#define QMCKL_NOT_PROVIDED ((qmckl_exit_code) 106)
|
||||
#define QMCKL_OUT_OF_BOUNDS ((qmckl_exit_code) 107)
|
||||
#define QMCKL_INVALID_EXIT_CODE ((qmckl_exit_code) 108)
|
||||
#define QMCKL_ALREADY_SET ((qmckl_exit_code) 108)
|
||||
#define QMCKL_INVALID_EXIT_CODE ((qmckl_exit_code) 109)
|
||||
#+end_src
|
||||
|
||||
#+begin_src f90 :comments org :tangle (eval fh_type) :exports none
|
||||
@ -196,7 +198,8 @@ return '\n'.join(result)
|
||||
integer(qmckl_exit_code), parameter :: QMCKL_DEALLOCATION_FAILED = 105
|
||||
integer(qmckl_exit_code), parameter :: QMCKL_NOT_PROVIDED = 106
|
||||
integer(qmckl_exit_code), parameter :: QMCKL_OUT_OF_BOUNDS = 107
|
||||
integer(qmckl_exit_code), parameter :: QMCKL_INVALID_EXIT_CODE = 108
|
||||
integer(qmckl_exit_code), parameter :: QMCKL_ALREADY_SET = 108
|
||||
integer(qmckl_exit_code), parameter :: QMCKL_INVALID_EXIT_CODE = 109
|
||||
#+end_src
|
||||
:end:
|
||||
|
||||
@ -239,7 +242,7 @@ for (text, code, message) in table:
|
||||
message = message.replace("'",'"')
|
||||
result += [ f"""case {text}:
|
||||
return {message};
|
||||
break;""" ]
|
||||
""" ]
|
||||
return '\n'.join(result)
|
||||
|
||||
#+end_src
|
||||
@ -248,88 +251,93 @@ return '\n'.join(result)
|
||||
#+begin_example
|
||||
case QMCKL_SUCCESS:
|
||||
return "Success";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_1:
|
||||
return "Invalid argument 1";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_2:
|
||||
return "Invalid argument 2";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_3:
|
||||
return "Invalid argument 3";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_4:
|
||||
return "Invalid argument 4";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_5:
|
||||
return "Invalid argument 5";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_6:
|
||||
return "Invalid argument 6";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_7:
|
||||
return "Invalid argument 7";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_8:
|
||||
return "Invalid argument 8";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_9:
|
||||
return "Invalid argument 9";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_10:
|
||||
return "Invalid argument 10";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_11:
|
||||
return "Invalid argument 11";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_12:
|
||||
return "Invalid argument 12";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_13:
|
||||
return "Invalid argument 13";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_14:
|
||||
return "Invalid argument 14";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_15:
|
||||
return "Invalid argument 15";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_16:
|
||||
return "Invalid argument 16";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_17:
|
||||
return "Invalid argument 17";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_18:
|
||||
return "Invalid argument 18";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_19:
|
||||
return "Invalid argument 19";
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_ARG_20:
|
||||
return "Invalid argument 20";
|
||||
break;
|
||||
|
||||
case QMCKL_FAILURE:
|
||||
return "Failure";
|
||||
break;
|
||||
|
||||
case QMCKL_ERRNO:
|
||||
return strerror(errno);
|
||||
break;
|
||||
|
||||
case QMCKL_INVALID_CONTEXT:
|
||||
return "Invalid context";
|
||||
break;
|
||||
|
||||
case QMCKL_ALLOCATION_FAILED:
|
||||
return "Allocation failed";
|
||||
break;
|
||||
|
||||
case QMCKL_DEALLOCATION_FAILED:
|
||||
return "De-allocation failed";
|
||||
break;
|
||||
|
||||
case QMCKL_NOT_PROVIDED:
|
||||
return "Not provided";
|
||||
break;
|
||||
|
||||
case QMCKL_OUT_OF_BOUNDS:
|
||||
return "Index out of bounds";
|
||||
|
||||
case QMCKL_ALREADY_SET:
|
||||
return "Already set";
|
||||
|
||||
case QMCKL_INVALID_EXIT_CODE:
|
||||
return "Invalid exit code";
|
||||
break;
|
||||
#+end_example
|
||||
|
||||
# Source
|
||||
@ -414,7 +422,7 @@ qmckl_set_error(qmckl_context context,
|
||||
|
||||
qmckl_lock(context);
|
||||
{
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL); /* Impossible because the context is valid. */
|
||||
|
||||
ctx->error.exit_code = exit_code;
|
||||
@ -460,7 +468,7 @@ qmckl_get_error(qmckl_context context,
|
||||
|
||||
qmckl_lock(context);
|
||||
{
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL); /* Impossible because the context is valid. */
|
||||
|
||||
/* Turn off annoying GCC warning */
|
||||
|
360
org/qmckl_examples.org
Normal file
360
org/qmckl_examples.org
Normal file
@ -0,0 +1,360 @@
|
||||
#+TITLE: Code examples
|
||||
#+SETUPFILE: ../tools/theme.setup
|
||||
#+INCLUDE: ../tools/lib.org
|
||||
|
||||
In this section, we present examples of usage of QMCkl.
|
||||
For simplicity, we assume that the wave function parameters are stored
|
||||
in a [[https://github.com/TREX-CoE/trexio][TREXIO]] file.
|
||||
|
||||
* Python
|
||||
** Check numerically that MOs are orthonormal
|
||||
|
||||
In this example, we will compute numerically the overlap
|
||||
between the molecular orbitals:
|
||||
|
||||
\[
|
||||
S_{ij} = \int \phi_i(\mathbf{r}) \phi_j(\mathbf{r})
|
||||
\text{d}\mathbf{r} \sim \sum_{k=1}^{N} \phi_i(\mathbf{r}_k)
|
||||
\phi_j(\mathbf{r}_k) \delta \mathbf{r}
|
||||
\]
|
||||
\[
|
||||
S_{ij} = \langle \phi_i | \phi_j \rangle
|
||||
\sim \sum_{k=1}^{N} \langle \phi_i | \mathbf{r}_k \rangle
|
||||
\langle \mathbf{r}_k | \phi_j \rangle
|
||||
\]
|
||||
|
||||
|
||||
#+begin_src python :exports code
|
||||
import numpy as np
|
||||
import qmckl
|
||||
#+end_src
|
||||
|
||||
#+RESULTS:
|
||||
|
||||
First, we create a context for the QMCkl calculation, and load the
|
||||
wave function stored in =h2o_5z.h5= inside it. It is a Hartree-Fock
|
||||
determinant for the water molecule in the cc-pV5Z basis set.
|
||||
|
||||
#+begin_src python :exports code
|
||||
trexio_filename = "..//share/qmckl/test_data/h2o_5z.h5"
|
||||
|
||||
context = qmckl.context_create()
|
||||
qmckl.trexio_read(context, trexio_filename)
|
||||
#+end_src
|
||||
|
||||
#+RESULTS:
|
||||
: None
|
||||
|
||||
We now define the grid points $\mathbf{r}_k$ as a regular grid around the
|
||||
molecule.
|
||||
|
||||
We fetch the nuclear coordinates from the context,
|
||||
|
||||
#+begin_src python :exports code
|
||||
nucl_num = qmckl.get_nucleus_num(context)
|
||||
|
||||
nucl_charge = qmckl.get_nucleus_charge(context, nucl_num)
|
||||
|
||||
nucl_coord = qmckl.get_nucleus_coord(context, 'N', nucl_num*3)
|
||||
nucl_coord = np.reshape(nucl_coord, (3, nucl_num))
|
||||
|
||||
for i in range(nucl_num):
|
||||
print("%d %+f %+f %+f"%(int(nucl_charge[i]),
|
||||
nucl_coord[i,0],
|
||||
nucl_coord[i,1],
|
||||
nucl_coord[i,2]) )
|
||||
#+end_src
|
||||
|
||||
#+begin_example
|
||||
8 +0.000000 +0.000000 +0.000000
|
||||
1 -1.430429 +0.000000 -1.107157
|
||||
1 +1.430429 +0.000000 -1.107157
|
||||
#+end_example
|
||||
|
||||
and compute the coordinates of the grid points:
|
||||
|
||||
#+begin_src python :exports code
|
||||
nx = ( 120, 120, 120 )
|
||||
shift = np.array([5.,5.,5.])
|
||||
point_num = nx[0] * nx[1] * nx[2]
|
||||
|
||||
rmin = np.array( list([ np.min(nucl_coord[:,a]) for a in range(3) ]) )
|
||||
rmax = np.array( list([ np.max(nucl_coord[:,a]) for a in range(3) ]) )
|
||||
|
||||
|
||||
linspace = [ None for i in range(3) ]
|
||||
step = [ None for i in range(3) ]
|
||||
for a in range(3):
|
||||
linspace[a], step[a] = np.linspace(rmin[a]-shift[a],
|
||||
rmax[a]+shift[a],
|
||||
num=nx[a],
|
||||
retstep=True)
|
||||
|
||||
dr = step[0] * step[1] * step[2]
|
||||
#+end_src
|
||||
|
||||
#+RESULTS:
|
||||
|
||||
Now the grid is ready, we can create the list of grid points
|
||||
$\mathbf{r}_k$ on which the MOs $\phi_i$ will be evaluated, and
|
||||
transfer them to the QMCkl context:
|
||||
|
||||
#+begin_src python :exports code
|
||||
point = []
|
||||
for x in linspace[0]:
|
||||
for y in linspace[1]:
|
||||
for z in linspace[2]:
|
||||
point += [ [x, y, z] ]
|
||||
|
||||
point = np.array(point)
|
||||
point_num = len(point)
|
||||
qmckl.set_point(context, 'N', point_num, np.reshape(point, (point_num*3)))
|
||||
#+end_src
|
||||
|
||||
#+RESULTS:
|
||||
: None
|
||||
|
||||
Then, we evaluate all the MOs at the grid points (and time the execution),
|
||||
and thus obtain the matrix $M_{ki} = \langle \mathbf{r}_k | \phi_i \rangle =
|
||||
\phi_i(\mathbf{r}_k)$.
|
||||
|
||||
#+begin_src python :exports code
|
||||
import time
|
||||
|
||||
mo_num = qmckl.get_mo_basis_mo_num(context)
|
||||
|
||||
before = time.time()
|
||||
mo_value = qmckl.get_mo_basis_mo_value(context, point_num*mo_num)
|
||||
after = time.time()
|
||||
|
||||
mo_value = np.reshape( mo_value, (point_num, mo_num) )
|
||||
|
||||
print("Number of MOs: ", mo_num)
|
||||
print("Number of grid points: ", point_num)
|
||||
print("Execution time : ", (after - before), "seconds")
|
||||
|
||||
#+end_src
|
||||
|
||||
#+begin_example
|
||||
Number of MOs: 201
|
||||
Number of grid points: 1728000
|
||||
Execution time : 3.511528968811035 seconds
|
||||
#+end_example
|
||||
|
||||
and finally we compute the overlap between all the MOs as
|
||||
$M^\dagger M$.
|
||||
|
||||
#+begin_src python :exports code
|
||||
overlap = mo_value.T @ mo_value * dr
|
||||
print (overlap)
|
||||
#+end_src
|
||||
|
||||
#+begin_example
|
||||
[[ 9.88693941e-01 2.34719693e-03 -1.50518232e-08 ... 3.12084178e-09
|
||||
-5.81064929e-10 3.70130091e-02]
|
||||
[ 2.34719693e-03 9.99509628e-01 3.18930040e-09 ... -2.46888958e-10
|
||||
-1.06064273e-09 -7.65567973e-03]
|
||||
[-1.50518232e-08 3.18930040e-09 9.99995073e-01 ... -5.84882580e-06
|
||||
-1.21598117e-06 4.59036468e-08]
|
||||
...
|
||||
[ 3.12084178e-09 -2.46888958e-10 -5.84882580e-06 ... 1.00019107e+00
|
||||
-2.03342837e-04 -1.36954855e-08]
|
||||
[-5.81064929e-10 -1.06064273e-09 -1.21598117e-06 ... -2.03342837e-04
|
||||
9.99262427e-01 1.18264754e-09]
|
||||
[ 3.70130091e-02 -7.65567973e-03 4.59036468e-08 ... -1.36954855e-08
|
||||
1.18264754e-09 8.97215950e-01]]
|
||||
#+end_example
|
||||
|
||||
* Fortran
|
||||
** Checking errors
|
||||
|
||||
All QMCkl functions return an error code. A convenient way to handle
|
||||
errors is to write an error-checking function that displays the
|
||||
error in text format and exits the program.
|
||||
|
||||
#+NAME: qmckl_check_error
|
||||
#+begin_src f90
|
||||
subroutine qmckl_check_error(rc, message)
|
||||
use qmckl
|
||||
implicit none
|
||||
integer(qmckl_exit_code), intent(in) :: rc
|
||||
character(len=*) , intent(in) :: message
|
||||
character(len=128) :: str_buffer
|
||||
if (rc /= QMCKL_SUCCESS) then
|
||||
print *, message
|
||||
call qmckl_string_of_error(rc, str_buffer)
|
||||
print *, str_buffer
|
||||
call exit(rc)
|
||||
end if
|
||||
end subroutine qmckl_check_error
|
||||
#+end_src
|
||||
|
||||
** Computing an atomic orbital on a grid
|
||||
:PROPERTIES:
|
||||
:header-args: :tangle ao_grid.f90
|
||||
:END:
|
||||
|
||||
The following program, in Fortran, computes the values of an atomic
|
||||
orbital on a regular 3-dimensional grid. The 100^3 grid points are
|
||||
automatically defined, such that the molecule fits in a box with 5
|
||||
atomic units in the borders.
|
||||
|
||||
This program uses the ~qmckl_check_error~ function defined above.
|
||||
|
||||
To use this program, run
|
||||
|
||||
#+begin_src bash :tangle no :exports code
|
||||
$ ao_grid <trexio_file> <AO_id> <point_num>
|
||||
#+end_src
|
||||
|
||||
|
||||
#+begin_src f90 :noweb yes
|
||||
<<qmckl_check_error>>
|
||||
|
||||
program ao_grid
|
||||
use qmckl
|
||||
implicit none
|
||||
|
||||
integer(qmckl_context) :: qmckl_ctx ! QMCkl context
|
||||
integer(qmckl_exit_code) :: rc ! Exit code of QMCkl functions
|
||||
|
||||
character(len=128) :: trexio_filename
|
||||
character(len=128) :: str_buffer
|
||||
integer :: ao_id
|
||||
integer :: point_num_x
|
||||
|
||||
integer(c_int64_t) :: nucl_num
|
||||
double precision, allocatable :: nucl_coord(:,:)
|
||||
|
||||
integer(c_int64_t) :: point_num
|
||||
integer(c_int64_t) :: ao_num
|
||||
integer(c_int64_t) :: ipoint, i, j, k
|
||||
double precision :: x, y, z, dr(3)
|
||||
double precision :: rmin(3), rmax(3)
|
||||
double precision, allocatable :: points(:,:)
|
||||
double precision, allocatable :: ao_vgl(:,:,:)
|
||||
#+end_src
|
||||
|
||||
Start by fetching the command-line arguments:
|
||||
|
||||
#+begin_src f90
|
||||
if (iargc() /= 3) then
|
||||
print *, 'Syntax: ao_grid <trexio_file> <AO_id> <point_num>'
|
||||
call exit(-1)
|
||||
end if
|
||||
call getarg(1, trexio_filename)
|
||||
call getarg(2, str_buffer)
|
||||
read(str_buffer, *) ao_id
|
||||
call getarg(3, str_buffer)
|
||||
read(str_buffer, *) point_num_x
|
||||
|
||||
if (point_num_x < 0 .or. point_num_x > 300) then
|
||||
print *, 'Error: 0 < point_num < 300'
|
||||
call exit(-1)
|
||||
end if
|
||||
#+end_src
|
||||
|
||||
Create the QMCkl context and initialize it with the wave function
|
||||
present in the TREXIO file:
|
||||
|
||||
#+begin_src f90
|
||||
qmckl_ctx = qmckl_context_create()
|
||||
rc = qmckl_trexio_read(qmckl_ctx, trexio_filename, 1_8*len(trim(trexio_filename)))
|
||||
call qmckl_check_error(rc, 'Read TREXIO')
|
||||
#+end_src
|
||||
|
||||
We need to check that ~ao_id~ is in the range, so we get the total
|
||||
number of AOs from QMCkl:
|
||||
|
||||
#+begin_src f90
|
||||
rc = qmckl_get_ao_basis_ao_num(qmckl_ctx, ao_num)
|
||||
call qmckl_check_error(rc, 'Getting ao_num')
|
||||
|
||||
if (ao_id < 0 .or. ao_id > ao_num) then
|
||||
print *, 'Error: 0 < ao_id < ', ao_num
|
||||
call exit(-1)
|
||||
end if
|
||||
#+end_src
|
||||
|
||||
Now we will compute the limits of the box in which the molecule fits.
|
||||
For that, we first need to ask QMCkl the coordinates of nuclei.
|
||||
|
||||
#+begin_src f90
|
||||
rc = qmckl_get_nucleus_num(qmckl_ctx, nucl_num)
|
||||
call qmckl_check_error(rc, 'Get nucleus num')
|
||||
|
||||
allocate( nucl_coord(3, nucl_num) )
|
||||
rc = qmckl_get_nucleus_coord(qmckl_ctx, 'N', nucl_coord, 3_8*nucl_num)
|
||||
call qmckl_check_error(rc, 'Get nucleus coord')
|
||||
#+end_src
|
||||
|
||||
We now compute the coordinates of opposite points of the box, and
|
||||
the distance between points along the 3 directions:
|
||||
|
||||
#+begin_src f90
|
||||
rmin(1) = minval( nucl_coord(1,:) ) - 5.d0
|
||||
rmin(2) = minval( nucl_coord(2,:) ) - 5.d0
|
||||
rmin(3) = minval( nucl_coord(3,:) ) - 5.d0
|
||||
|
||||
rmax(1) = maxval( nucl_coord(1,:) ) + 5.d0
|
||||
rmax(2) = maxval( nucl_coord(2,:) ) + 5.d0
|
||||
rmax(3) = maxval( nucl_coord(3,:) ) + 5.d0
|
||||
|
||||
dr(1:3) = (rmax(1:3) - rmin(1:3)) / dble(point_num_x-1)
|
||||
#+end_src
|
||||
|
||||
We now produce the list of point coordinates where the AO will be
|
||||
evaluated:
|
||||
|
||||
#+begin_src f90
|
||||
point_num = point_num_x**3
|
||||
allocate( points(point_num, 3) )
|
||||
ipoint=0
|
||||
z = rmin(3)
|
||||
do k=1,point_num_x
|
||||
y = rmin(2)
|
||||
do j=1,point_num_x
|
||||
x = rmin(1)
|
||||
do i=1,point_num_x
|
||||
ipoint = ipoint+1
|
||||
points(ipoint,1) = x
|
||||
points(ipoint,2) = y
|
||||
points(ipoint,3) = z
|
||||
x = x + dr(1)
|
||||
end do
|
||||
y = y + dr(2)
|
||||
end do
|
||||
z = z + dr(3)
|
||||
end do
|
||||
#+end_src
|
||||
|
||||
We give the points to QMCkl:
|
||||
|
||||
#+begin_src f90
|
||||
rc = qmckl_set_point(qmckl_ctx, 'T', point_num, points, size(points)*1_8 )
|
||||
call qmckl_check_error(rc, 'Setting points')
|
||||
#+end_src
|
||||
|
||||
We allocate the space required to retrieve the values, gradients and
|
||||
Laplacian of all AOs, and ask to retrieve the values of the
|
||||
AOs computed at the point positions.
|
||||
|
||||
#+begin_src f90
|
||||
allocate( ao_vgl(ao_num, 5, point_num) )
|
||||
rc = qmckl_get_ao_basis_ao_vgl(qmckl_ctx, ao_vgl, ao_num*5_8*point_num)
|
||||
call qmckl_check_error(rc, 'Setting points')
|
||||
#+end_src
|
||||
|
||||
We finally print the value and Laplacian of the AO:
|
||||
|
||||
#+begin_src f90
|
||||
do ipoint=1, point_num
|
||||
print '(3(F10.6,X),2(E20.10,X))', points(ipoint, 1:3), ao_vgl(ao_id,1,ipoint), ao_vgl(ao_id,5,ipoint)
|
||||
end do
|
||||
#+end_src
|
||||
|
||||
#+begin_src f90
|
||||
deallocate( nucl_coord, points, ao_vgl )
|
||||
end program ao_grid
|
||||
#+end_src
|
||||
|
File diff suppressed because it is too large
Load Diff
@ -226,7 +226,7 @@ qmckl_exit_code qmckl_get_kinetic_energy(qmckl_context context, double * const k
|
||||
rc = qmckl_provide_kinetic_energy(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->electron.walk_num * sizeof(double);
|
||||
@ -250,7 +250,7 @@ qmckl_exit_code qmckl_provide_kinetic_energy(qmckl_context context) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if(!(ctx->nucleus.provided)) {
|
||||
@ -549,37 +549,28 @@ end function qmckl_compute_kinetic_energy_f
|
||||
|
||||
*** Test
|
||||
#+begin_src c :tangle (eval c_test) :exports none
|
||||
|
||||
#define walk_num chbrclf_walk_num
|
||||
#define elec_num chbrclf_elec_num
|
||||
#define shell_num chbrclf_shell_num
|
||||
#define ao_num chbrclf_ao_num
|
||||
|
||||
int64_t elec_up_num = chbrclf_elec_up_num;
|
||||
int64_t elec_dn_num = chbrclf_elec_dn_num;
|
||||
double* elec_coord = &(chbrclf_elec_coord[0][0][0]);
|
||||
const int64_t nucl_num = chbrclf_nucl_num;
|
||||
const double* nucl_charge = chbrclf_charge;
|
||||
const double* nucl_coord = &(chbrclf_nucl_coord[0][0]);
|
||||
|
||||
rc = qmckl_set_electron_num (context, elec_up_num, elec_dn_num);
|
||||
rc = qmckl_set_electron_num (context, chbrclf_elec_up_num, chbrclf_elec_dn_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_electron_walk_num (context, walk_num);
|
||||
rc = qmckl_set_electron_walk_num (context, chbrclf_walk_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
assert(qmckl_electron_provided(context));
|
||||
|
||||
rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*elec_num*3);
|
||||
rc = qmckl_set_electron_coord (context, 'N', elec_coord, chbrclf_walk_num*chbrclf_elec_num*3);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_num (context, nucl_num);
|
||||
rc = qmckl_set_nucleus_num (context, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0]), nucl_num*3);
|
||||
rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0]), chbrclf_nucl_num*3);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_nucleus_charge(context, nucl_charge, nucl_num);
|
||||
rc = qmckl_set_nucleus_charge(context, nucl_charge, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
assert(qmckl_nucleus_provided(context));
|
||||
@ -611,11 +602,11 @@ rc = qmckl_set_ao_basis_prim_num (context, chbrclf_prim_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_ao_basis_provided(context));
|
||||
|
||||
rc = qmckl_set_ao_basis_nucleus_index (context, nucleus_index, nucl_num);
|
||||
rc = qmckl_set_ao_basis_nucleus_index (context, nucleus_index, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_ao_basis_provided(context));
|
||||
|
||||
rc = qmckl_set_ao_basis_nucleus_shell_num (context, nucleus_shell_num, nucl_num);
|
||||
rc = qmckl_set_ao_basis_nucleus_shell_num (context, nucleus_shell_num, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_ao_basis_provided(context));
|
||||
|
||||
@ -655,10 +646,10 @@ assert(rc == QMCKL_SUCCESS);
|
||||
assert(qmckl_ao_basis_provided(context));
|
||||
|
||||
|
||||
double ao_vgl[walk_num*elec_num][5][chbrclf_ao_num];
|
||||
double ao_vgl[chbrclf_walk_num*chbrclf_elec_num][5][chbrclf_ao_num];
|
||||
|
||||
rc = qmckl_get_ao_basis_ao_vgl(context, &(ao_vgl[0][0][0]),
|
||||
(int64_t) 5*walk_num*elec_num*chbrclf_ao_num);
|
||||
(int64_t) 5*chbrclf_walk_num*chbrclf_elec_num*chbrclf_ao_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
/* Set up MO data */
|
||||
@ -673,31 +664,31 @@ assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
assert(qmckl_mo_basis_provided(context));
|
||||
|
||||
double mo_vgl[walk_num*elec_num][5][chbrclf_mo_num];
|
||||
rc = qmckl_get_mo_basis_mo_vgl(context, &(mo_vgl[0][0][0]), 5*walk_num*elec_num*chbrclf_mo_num);
|
||||
double mo_vgl[chbrclf_walk_num*chbrclf_elec_num][5][chbrclf_mo_num];
|
||||
rc = qmckl_get_mo_basis_mo_vgl(context, &(mo_vgl[0][0][0]), 5*chbrclf_walk_num*chbrclf_elec_num*chbrclf_mo_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
/* Set up determinant data */
|
||||
|
||||
const int64_t det_num_alpha = 1;
|
||||
const int64_t det_num_beta = 1;
|
||||
int64_t mo_index_alpha[det_num_alpha][walk_num][elec_up_num];
|
||||
int64_t mo_index_beta[det_num_alpha][walk_num][elec_dn_num];
|
||||
#define det_num_alpha 1
|
||||
#define det_num_beta 1
|
||||
int64_t mo_index_alpha[det_num_alpha][chbrclf_walk_num][chbrclf_elec_up_num];
|
||||
int64_t mo_index_beta[det_num_alpha][chbrclf_walk_num][chbrclf_elec_dn_num];
|
||||
|
||||
int i, j, k;
|
||||
for(k = 0; k < det_num_alpha; ++k)
|
||||
for(i = 0; i < walk_num; ++i)
|
||||
for(j = 0; j < elec_up_num; ++j)
|
||||
for(i = 0; i < chbrclf_walk_num; ++i)
|
||||
for(j = 0; j < chbrclf_elec_up_num; ++j)
|
||||
mo_index_alpha[k][i][j] = j + 1;
|
||||
for(k = 0; k < det_num_beta; ++k)
|
||||
for(i = 0; i < walk_num; ++i)
|
||||
for(j = 0; j < elec_up_num; ++j)
|
||||
for(i = 0; i < chbrclf_walk_num; ++i)
|
||||
for(j = 0; j < chbrclf_elec_up_num; ++j)
|
||||
mo_index_beta[k][i][j] = j + 1;
|
||||
|
||||
rc = qmckl_set_determinant_type (context, typ);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_determinant_walk_num (context, walk_num);
|
||||
rc = qmckl_set_determinant_walk_num (context, chbrclf_walk_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_set_determinant_det_num_alpha (context, det_num_alpha);
|
||||
@ -714,8 +705,8 @@ assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
// Get alpha determinant
|
||||
|
||||
double det_vgl_alpha[det_num_alpha][walk_num][5][elec_up_num][elec_up_num];
|
||||
double det_vgl_beta[det_num_beta][walk_num][5][elec_dn_num][elec_dn_num];
|
||||
double det_vgl_alpha[det_num_alpha][chbrclf_walk_num][5][chbrclf_elec_up_num][chbrclf_elec_up_num];
|
||||
double det_vgl_beta[det_num_beta][chbrclf_walk_num][5][chbrclf_elec_dn_num][chbrclf_elec_dn_num];
|
||||
|
||||
rc = qmckl_get_det_vgl_alpha(context, &(det_vgl_alpha[0][0][0][0][0]));
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
@ -725,8 +716,8 @@ assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
// Get adjoint of the slater-determinant
|
||||
|
||||
double det_inv_matrix_alpha[det_num_alpha][walk_num][elec_up_num][elec_up_num];
|
||||
double det_inv_matrix_beta[det_num_beta][walk_num][elec_dn_num][elec_dn_num];
|
||||
double det_inv_matrix_alpha[det_num_alpha][chbrclf_walk_num][chbrclf_elec_up_num][chbrclf_elec_up_num];
|
||||
double det_inv_matrix_beta[det_num_beta][chbrclf_walk_num][chbrclf_elec_dn_num][chbrclf_elec_dn_num];
|
||||
|
||||
rc = qmckl_get_det_inv_matrix_alpha(context, &(det_inv_matrix_alpha[0][0][0][0]));
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
@ -736,7 +727,7 @@ assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
// Calculate the Kinetic energy
|
||||
|
||||
double kinetic_energy[walk_num];
|
||||
double kinetic_energy[chbrclf_walk_num];
|
||||
|
||||
rc = qmckl_get_kinetic_energy(context, &(kinetic_energy[0]));
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
@ -799,7 +790,7 @@ qmckl_exit_code qmckl_get_potential_energy(qmckl_context context, double * const
|
||||
rc = qmckl_provide_potential_energy(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->electron.walk_num * sizeof(double);
|
||||
@ -822,7 +813,7 @@ qmckl_exit_code qmckl_provide_potential_energy(qmckl_context context) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
qmckl_exit_code rc;
|
||||
@ -1034,7 +1025,7 @@ end function qmckl_compute_potential_energy_f
|
||||
#+begin_src c :tangle (eval c_test) :exports none
|
||||
// Calculate the Potential energy
|
||||
|
||||
double potential_energy[walk_num];
|
||||
double potential_energy[chbrclf_walk_num];
|
||||
|
||||
rc = qmckl_get_potential_energy(context, &(potential_energy[0]));
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
@ -1083,7 +1074,7 @@ qmckl_exit_code qmckl_get_local_energy(qmckl_context context, double * const loc
|
||||
rc = qmckl_provide_local_energy(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
const int64_t sze = ctx->electron.walk_num;
|
||||
@ -1112,7 +1103,7 @@ qmckl_exit_code qmckl_provide_local_energy(qmckl_context context) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if(!(ctx->nucleus.provided)) {
|
||||
@ -1296,9 +1287,9 @@ end function qmckl_compute_local_energy_f
|
||||
#+begin_src c :tangle (eval c_test) :exports none
|
||||
// Calculate the Local energy
|
||||
|
||||
double local_energy[walk_num];
|
||||
double local_energy[chbrclf_walk_num];
|
||||
|
||||
rc = qmckl_get_local_energy(context, &(local_energy[0]), walk_num);
|
||||
rc = qmckl_get_local_energy(context, &(local_energy[0]), chbrclf_walk_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
#+end_src
|
||||
@ -1345,7 +1336,7 @@ qmckl_exit_code qmckl_get_drift_vector(qmckl_context context, double * const dri
|
||||
rc = qmckl_provide_drift_vector(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
size_t sze = ctx->electron.walk_num * ctx->electron.num * 3 * sizeof(double);
|
||||
@ -1368,7 +1359,7 @@ qmckl_exit_code qmckl_provide_drift_vector(qmckl_context context) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if(!(ctx->nucleus.provided)) {
|
||||
@ -1651,7 +1642,7 @@ end function qmckl_compute_drift_vector_f
|
||||
#+begin_src c :tangle (eval c_test) :exports none
|
||||
// Calculate the Drift vector
|
||||
|
||||
double drift_vector[walk_num][elec_num][3];
|
||||
double drift_vector[chbrclf_walk_num][chbrclf_elec_num][3];
|
||||
|
||||
rc = qmckl_get_drift_vector(context, &(drift_vector[0][0][0]));
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
@ -116,10 +116,14 @@ void* qmckl_malloc(qmckl_context context, const qmckl_memory_info_struct info) {
|
||||
|
||||
assert (qmckl_context_check(context) != QMCKL_NULL_CONTEXT);
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
|
||||
/* Allocate memory and zero it */
|
||||
#ifdef HAVE_HPC
|
||||
void * pointer = aligned_alloc(64, ((info.size+64) >> 6) << 6 );
|
||||
#else
|
||||
void * pointer = malloc(info.size);
|
||||
#endif
|
||||
if (pointer == NULL) {
|
||||
return NULL;
|
||||
}
|
||||
@ -217,7 +221,7 @@ qmckl_exit_code qmckl_free(qmckl_context context, void * const ptr) {
|
||||
"NULL pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
|
||||
qmckl_lock(context);
|
||||
{
|
||||
|
609
org/qmckl_mo.org
609
org/qmckl_mo.org
@ -92,10 +92,12 @@ int main() {
|
||||
|
||||
Computed data:
|
||||
|
||||
|---------------+--------------------------+-------------------------------------------------------------------------------------|
|
||||
|-----------------+--------------------------+-------------------------------------------------------------------------------------|
|
||||
| ~mo_value~ | ~[point_num][mo_num]~ | Value of the MOs at point positions |
|
||||
| ~mo_value_date~ | ~uint64_t~ | Late modification date of the value of the MOs at point positions |
|
||||
| ~mo_vgl~ | ~[point_num][5][mo_num]~ | Value, gradients, Laplacian of the MOs at point positions |
|
||||
| ~mo_vgl_date~ | ~uint64_t~ | Late modification date of Value, gradients, Laplacian of the MOs at point positions |
|
||||
|---------------+--------------------------+-------------------------------------------------------------------------------------|
|
||||
|-----------------+--------------------------+-------------------------------------------------------------------------------------|
|
||||
|
||||
** Data structure
|
||||
|
||||
@ -106,7 +108,9 @@ typedef struct qmckl_mo_basis_struct {
|
||||
double * restrict coefficient_t;
|
||||
|
||||
double * restrict mo_vgl;
|
||||
double * restrict mo_value;
|
||||
uint64_t mo_vgl_date;
|
||||
uint64_t mo_value_date;
|
||||
|
||||
int32_t uninitialized;
|
||||
bool provided;
|
||||
@ -131,7 +135,7 @@ qmckl_exit_code qmckl_init_mo_basis(qmckl_context context) {
|
||||
return false;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
ctx->mo_basis.uninitialized = (1 << 2) - 1;
|
||||
@ -158,10 +162,9 @@ qmckl_get_mo_basis_mo_num (const qmckl_context context,
|
||||
QMCKL_INVALID_CONTEXT,
|
||||
"qmckl_get_mo_basis_mo_num",
|
||||
NULL);
|
||||
return (int64_t) 0;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1;
|
||||
@ -200,7 +203,7 @@ qmckl_get_mo_basis_coefficient (const qmckl_context context,
|
||||
NULL);
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 1;
|
||||
@ -248,7 +251,7 @@ bool qmckl_mo_basis_provided(const qmckl_context context) {
|
||||
return false;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
return ctx->mo_basis.provided;
|
||||
@ -257,10 +260,9 @@ bool qmckl_mo_basis_provided(const qmckl_context context) {
|
||||
|
||||
#+end_src
|
||||
|
||||
|
||||
*** Fortran interfaces
|
||||
|
||||
#+begin_src f90 :tangle (eval fh_func) :comments org
|
||||
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
||||
interface
|
||||
integer(c_int32_t) function qmckl_get_mo_basis_mo_num (context, &
|
||||
mo_num) bind(C)
|
||||
@ -280,7 +282,7 @@ interface
|
||||
implicit none
|
||||
integer (c_int64_t) , intent(in) , value :: context
|
||||
double precision, intent(out) :: coefficient(*)
|
||||
integer (c_int64_t) , intent(in) , value :: size_max
|
||||
integer (c_int64_t) , intent(in), value :: size_max
|
||||
end function qmckl_get_mo_basis_coefficient
|
||||
end interface
|
||||
|
||||
@ -302,7 +304,14 @@ if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
|
||||
if (mask != 0 && !(ctx->mo_basis.uninitialized & mask)) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_ALREADY_SET,
|
||||
"qmckl_set_mo_*",
|
||||
NULL);
|
||||
}
|
||||
#+end_src
|
||||
|
||||
#+NAME:post
|
||||
@ -318,6 +327,9 @@ return QMCKL_SUCCESS;
|
||||
|
||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
||||
qmckl_exit_code qmckl_set_mo_basis_mo_num(qmckl_context context, const int64_t mo_num) {
|
||||
|
||||
int32_t mask = 1 ;
|
||||
|
||||
<<pre>>
|
||||
|
||||
if (mo_num <= 0) {
|
||||
@ -327,17 +339,17 @@ qmckl_exit_code qmckl_set_mo_basis_mo_num(qmckl_context context, const int64_t m
|
||||
"mo_num <= 0");
|
||||
}
|
||||
|
||||
int32_t mask = 1 ;
|
||||
ctx->mo_basis.mo_num = mo_num;
|
||||
|
||||
<<post>>
|
||||
}
|
||||
|
||||
qmckl_exit_code qmckl_set_mo_basis_coefficient(qmckl_context context, const double* coefficient) {
|
||||
<<pre>>
|
||||
|
||||
int32_t mask = 1 << 1;
|
||||
|
||||
<<pre>>
|
||||
|
||||
if (ctx->mo_basis.coefficient != NULL) {
|
||||
qmckl_exit_code rc = qmckl_free(context, ctx->mo_basis.coefficient);
|
||||
if (rc != QMCKL_SUCCESS) {
|
||||
@ -383,7 +395,7 @@ qmckl_exit_code qmckl_finalize_mo_basis(qmckl_context context) {
|
||||
NULL);
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
|
||||
@ -421,7 +433,464 @@ qmckl_exit_code qmckl_finalize_mo_basis(qmckl_context context) {
|
||||
|
||||
* Computation
|
||||
|
||||
** Computation of MOs
|
||||
** Computation of MOs: values only
|
||||
|
||||
*** Get
|
||||
|
||||
#+begin_src c :comments org :tangle (eval h_func) :noweb yes
|
||||
qmckl_exit_code
|
||||
qmckl_get_mo_basis_mo_value(qmckl_context context,
|
||||
double* const mo_value,
|
||||
const int64_t size_max);
|
||||
#+end_src
|
||||
|
||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
||||
qmckl_exit_code
|
||||
qmckl_get_mo_basis_mo_value(qmckl_context context,
|
||||
double* const mo_value,
|
||||
const int64_t size_max)
|
||||
{
|
||||
|
||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_exit_code rc;
|
||||
|
||||
rc = qmckl_provide_ao_basis_ao_value(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
rc = qmckl_provide_mo_value(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
const int64_t sze = ctx->point.num * ctx->mo_basis.mo_num;
|
||||
if (size_max < sze) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_INVALID_ARG_3,
|
||||
"qmckl_get_mo_basis_mo_value",
|
||||
"input array too small");
|
||||
}
|
||||
memcpy(mo_value, ctx->mo_basis.mo_value, sze * sizeof(double));
|
||||
|
||||
return QMCKL_SUCCESS;
|
||||
}
|
||||
#+end_src
|
||||
|
||||
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
||||
interface
|
||||
integer(c_int32_t) function qmckl_get_mo_basis_mo_value (context, &
|
||||
mo_value, size_max) bind(C)
|
||||
use, intrinsic :: iso_c_binding
|
||||
import
|
||||
implicit none
|
||||
|
||||
integer (c_int64_t) , intent(in) , value :: context
|
||||
double precision, intent(out) :: mo_value(*)
|
||||
integer (c_int64_t) , intent(in) , value :: size_max
|
||||
end function qmckl_get_mo_basis_mo_value
|
||||
end interface
|
||||
#+end_src
|
||||
|
||||
Uses the given array to compute the values.
|
||||
|
||||
#+begin_src c :comments org :tangle (eval h_func) :noweb yes
|
||||
qmckl_exit_code
|
||||
qmckl_get_mo_basis_mo_value_inplace (qmckl_context context,
|
||||
double* const mo_value,
|
||||
const int64_t size_max);
|
||||
#+end_src
|
||||
|
||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
||||
qmckl_exit_code
|
||||
qmckl_get_mo_basis_mo_value_inplace (qmckl_context context,
|
||||
double* const mo_value,
|
||||
const int64_t size_max)
|
||||
{
|
||||
|
||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_INVALID_CONTEXT,
|
||||
"qmckl_get_mo_basis_mo_value",
|
||||
NULL);
|
||||
}
|
||||
|
||||
qmckl_exit_code rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
const int64_t sze = ctx->mo_basis.mo_num * ctx->point.num;
|
||||
if (size_max < sze) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_INVALID_ARG_3,
|
||||
"qmckl_get_mo_basis_mo_value",
|
||||
"input array too small");
|
||||
}
|
||||
|
||||
rc = qmckl_context_touch(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
double* old_array = ctx->mo_basis.mo_value;
|
||||
|
||||
ctx->mo_basis.mo_value = mo_value;
|
||||
|
||||
rc = qmckl_provide_mo_value(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
ctx->mo_basis.mo_value = old_array;
|
||||
|
||||
return QMCKL_SUCCESS;
|
||||
}
|
||||
#+end_src
|
||||
|
||||
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
||||
interface
|
||||
integer(c_int32_t) function qmckl_get_mo_basis_mo_value_inplace (context, &
|
||||
mo_value, size_max) bind(C)
|
||||
use, intrinsic :: iso_c_binding
|
||||
import
|
||||
implicit none
|
||||
integer (c_int64_t) , intent(in) , value :: context
|
||||
double precision, intent(out) :: mo_value(*)
|
||||
integer (c_int64_t) , intent(in) , value :: size_max
|
||||
end function qmckl_get_mo_basis_mo_value_inplace
|
||||
end interface
|
||||
#+end_src
|
||||
|
||||
*** Provide
|
||||
|
||||
#+begin_src c :comments org :tangle (eval h_private_func) :noweb yes :exports none
|
||||
qmckl_exit_code qmckl_provide_mo_value(qmckl_context context);
|
||||
#+end_src
|
||||
|
||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
||||
qmckl_exit_code qmckl_provide_mo_value(qmckl_context context)
|
||||
{
|
||||
|
||||
qmckl_exit_code rc;
|
||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!ctx->ao_basis.provided) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_NOT_PROVIDED,
|
||||
"qmckl_ao_basis",
|
||||
NULL);
|
||||
}
|
||||
|
||||
rc = qmckl_provide_ao_basis_ao_value(context);
|
||||
if (rc != QMCKL_SUCCESS) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_NOT_PROVIDED,
|
||||
"qmckl_ao_value",
|
||||
NULL);
|
||||
}
|
||||
|
||||
if (!ctx->mo_basis.provided) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_NOT_PROVIDED,
|
||||
"qmckl_mo_basis",
|
||||
NULL);
|
||||
}
|
||||
|
||||
/* Compute if necessary */
|
||||
if (ctx->point.date > ctx->mo_basis.mo_value_date) {
|
||||
|
||||
/* Allocate array */
|
||||
if (ctx->mo_basis.mo_value == NULL) {
|
||||
|
||||
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
|
||||
mem_info.size = ctx->point.num * ctx->mo_basis.mo_num * sizeof(double);
|
||||
double* mo_value = (double*) qmckl_malloc(context, mem_info);
|
||||
|
||||
if (mo_value == NULL) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_ALLOCATION_FAILED,
|
||||
"qmckl_mo_basis_mo_value",
|
||||
NULL);
|
||||
}
|
||||
ctx->mo_basis.mo_value = mo_value;
|
||||
}
|
||||
|
||||
if (ctx->mo_basis.mo_vgl_date == ctx->point.date) {
|
||||
|
||||
// mo_vgl has been computed at this step: Just copy the data.
|
||||
|
||||
double * v = &(ctx->mo_basis.mo_value[0]);
|
||||
double * vgl = &(ctx->mo_basis.mo_vgl[0]);
|
||||
for (int i=0 ; i<ctx->point.num ; ++i) {
|
||||
for (int k=0 ; k<ctx->mo_basis.mo_num ; ++k) {
|
||||
v[k] = vgl[k];
|
||||
}
|
||||
v += ctx->mo_basis.mo_num;
|
||||
vgl += ctx->mo_basis.mo_num * 5;
|
||||
}
|
||||
|
||||
} else {
|
||||
|
||||
rc = qmckl_compute_mo_basis_mo_value(context,
|
||||
ctx->ao_basis.ao_num,
|
||||
ctx->mo_basis.mo_num,
|
||||
ctx->point.num,
|
||||
ctx->mo_basis.coefficient_t,
|
||||
ctx->ao_basis.ao_value,
|
||||
ctx->mo_basis.mo_value);
|
||||
|
||||
if (rc != QMCKL_SUCCESS) {
|
||||
return rc;
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
ctx->mo_basis.mo_value_date = ctx->date;
|
||||
}
|
||||
|
||||
return QMCKL_SUCCESS;
|
||||
}
|
||||
#+end_src
|
||||
|
||||
*** Compute
|
||||
:PROPERTIES:
|
||||
:Name: qmckl_compute_mo_basis_mo_value
|
||||
:CRetType: qmckl_exit_code
|
||||
:FRetType: qmckl_exit_code
|
||||
:END:
|
||||
|
||||
#+NAME: qmckl_mo_basis_mo_value_args
|
||||
| Variable | Type | In/Out | Description |
|
||||
|-----------------+-----------------------------+--------+-------------------------------------------------|
|
||||
| ~context~ | ~qmckl_context~ | in | Global state |
|
||||
| ~ao_num~ | ~int64_t~ | in | Number of AOs |
|
||||
| ~mo_num~ | ~int64_t~ | in | Number of MOs |
|
||||
| ~point_num~ | ~int64_t~ | in | Number of points |
|
||||
| ~coefficient_t~ | ~double[mo_num][ao_num]~ | in | Transpose of the AO to MO transformation matrix |
|
||||
| ~ao_value~ | ~double[point_num][ao_num]~ | in | Value of the AOs |
|
||||
| ~mo_value~ | ~double[point_num][mo_num]~ | out | Value of the MOs |
|
||||
|
||||
|
||||
The matrix of AO values is very sparse, so we use a sparse-dense
|
||||
matrix multiplication instead of a dgemm, as exposed in
|
||||
https://dx.doi.org/10.1007/978-3-642-38718-0_14.
|
||||
|
||||
|
||||
|
||||
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
|
||||
integer function qmckl_compute_mo_basis_mo_value_doc_f(context, &
|
||||
ao_num, mo_num, point_num, &
|
||||
coefficient_t, ao_value, mo_value) &
|
||||
result(info)
|
||||
use qmckl
|
||||
implicit none
|
||||
integer(qmckl_context), intent(in) :: context
|
||||
integer*8 , intent(in) :: ao_num, mo_num
|
||||
integer*8 , intent(in) :: point_num
|
||||
double precision , intent(in) :: ao_value(ao_num,point_num)
|
||||
double precision , intent(in) :: coefficient_t(mo_num,ao_num)
|
||||
double precision , intent(out) :: mo_value(mo_num,point_num)
|
||||
integer*8 :: i,j,k
|
||||
double precision :: c1, c2, c3, c4, c5
|
||||
|
||||
integer*8 :: LDA, LDB, LDC
|
||||
|
||||
info = QMCKL_SUCCESS
|
||||
if (.True.) then ! fast algorithm
|
||||
do j=1,point_num
|
||||
mo_value(:,j) = 0.d0
|
||||
do k=1,ao_num
|
||||
if (ao_value(k,j) /= 0.d0) then
|
||||
c1 = ao_value(k,j)
|
||||
do i=1,mo_num
|
||||
mo_value(i,j) = mo_value(i,j) + coefficient_t(i,k) * c1
|
||||
end do
|
||||
end if
|
||||
end do
|
||||
end do
|
||||
|
||||
else ! dgemm
|
||||
|
||||
LDA = size(coefficient_t,1)
|
||||
LDB = size(ao_value,1)
|
||||
LDC = size(mo_value,1)
|
||||
|
||||
info = qmckl_dgemm(context,'N', 'N', mo_num, point_num, ao_num, 1.d0, &
|
||||
coefficient_t, LDA, ao_value, LDB, &
|
||||
0.d0, mo_value, LDC)
|
||||
|
||||
end if
|
||||
|
||||
end function qmckl_compute_mo_basis_mo_value_doc_f
|
||||
#+end_src
|
||||
|
||||
#+CALL: generate_c_header(table=qmckl_mo_basis_mo_value_args,rettyp=get_value("CRetType"),fname="qmckl_compute_mo_basis_mo_value"))
|
||||
|
||||
#+RESULTS:
|
||||
#+begin_src c :tangle (eval h_func) :comments org
|
||||
qmckl_exit_code qmckl_compute_mo_basis_mo_value (
|
||||
const qmckl_context context,
|
||||
const int64_t ao_num,
|
||||
const int64_t mo_num,
|
||||
const int64_t point_num,
|
||||
const double* coefficient_t,
|
||||
const double* ao_value,
|
||||
double* const mo_value );
|
||||
#+end_src
|
||||
|
||||
#+CALL: generate_c_header(table=qmckl_mo_basis_mo_value_args,rettyp=get_value("CRetType"),fname="qmckl_compute_mo_basis_mo_value_doc"))
|
||||
|
||||
#+RESULTS:
|
||||
#+begin_src c :tangle (eval h_func) :comments org
|
||||
qmckl_exit_code qmckl_compute_mo_basis_mo_value_doc (
|
||||
const qmckl_context context,
|
||||
const int64_t ao_num,
|
||||
const int64_t mo_num,
|
||||
const int64_t point_num,
|
||||
const double* coefficient_t,
|
||||
const double* ao_value,
|
||||
double* const mo_value );
|
||||
#+end_src
|
||||
|
||||
#+CALL: generate_c_interface(table=qmckl_mo_basis_mo_value_args,rettyp=get_value("CRetType"),fname="qmckl_compute_mo_basis_mo_value_doc"))
|
||||
|
||||
#+RESULTS:
|
||||
#+begin_src f90 :tangle (eval f) :comments org :exports none
|
||||
integer(c_int32_t) function qmckl_compute_mo_basis_mo_value_doc &
|
||||
(context, ao_num, mo_num, point_num, coefficient_t, ao_value, mo_value) &
|
||||
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 :: ao_num
|
||||
integer (c_int64_t) , intent(in) , value :: mo_num
|
||||
integer (c_int64_t) , intent(in) , value :: point_num
|
||||
real (c_double ) , intent(in) :: coefficient_t(ao_num,mo_num)
|
||||
real (c_double ) , intent(in) :: ao_value(ao_num,point_num)
|
||||
real (c_double ) , intent(out) :: mo_value(mo_num,point_num)
|
||||
|
||||
integer(c_int32_t), external :: qmckl_compute_mo_basis_mo_value_doc_f
|
||||
info = qmckl_compute_mo_basis_mo_value_doc_f &
|
||||
(context, ao_num, mo_num, point_num, coefficient_t, ao_value, mo_value)
|
||||
|
||||
end function qmckl_compute_mo_basis_mo_value_doc
|
||||
#+end_src
|
||||
|
||||
#+begin_src c :tangle (eval c) :comments org
|
||||
qmckl_exit_code
|
||||
qmckl_compute_mo_basis_mo_value (const qmckl_context context,
|
||||
const int64_t ao_num,
|
||||
const int64_t mo_num,
|
||||
const int64_t point_num,
|
||||
const double* coefficient_t,
|
||||
const double* ao_value,
|
||||
double* const mo_value )
|
||||
{
|
||||
#ifdef HAVE_HPC
|
||||
return qmckl_compute_mo_basis_mo_value_hpc (context, ao_num, mo_num, point_num, coefficient_t, ao_value, mo_value);
|
||||
#else
|
||||
return qmckl_compute_mo_basis_mo_value_doc (context, ao_num, mo_num, point_num, coefficient_t, ao_value, mo_value);
|
||||
#endif
|
||||
}
|
||||
#+end_src
|
||||
|
||||
*** HPC version
|
||||
|
||||
|
||||
#+begin_src c :tangle (eval h_func) :comments org
|
||||
#ifdef HAVE_HPC
|
||||
qmckl_exit_code
|
||||
qmckl_compute_mo_basis_mo_value_hpc (const qmckl_context context,
|
||||
const int64_t ao_num,
|
||||
const int64_t mo_num,
|
||||
const int64_t point_num,
|
||||
const double* coefficient_t,
|
||||
const double* ao_value,
|
||||
double* const mo_value );
|
||||
#endif
|
||||
#+end_src
|
||||
|
||||
#+begin_src c :tangle (eval c) :comments org
|
||||
#ifdef HAVE_HPC
|
||||
qmckl_exit_code
|
||||
qmckl_compute_mo_basis_mo_value_hpc (const qmckl_context context,
|
||||
const int64_t ao_num,
|
||||
const int64_t mo_num,
|
||||
const int64_t point_num,
|
||||
const double* restrict coefficient_t,
|
||||
const double* restrict ao_value,
|
||||
double* restrict const mo_value )
|
||||
{
|
||||
assert (context != QMCKL_NULL_CONTEXT);
|
||||
|
||||
#ifdef HAVE_OPENMP
|
||||
#pragma omp parallel for
|
||||
#endif
|
||||
for (int64_t ipoint=0 ; ipoint < point_num ; ++ipoint) {
|
||||
double* restrict const vgl1 = &(mo_value[ipoint*mo_num]);
|
||||
const double* restrict avgl1 = &(ao_value[ipoint*ao_num]);
|
||||
|
||||
for (int64_t i=0 ; i<mo_num ; ++i) {
|
||||
vgl1[i] = 0.;
|
||||
}
|
||||
|
||||
int64_t nidx=0;
|
||||
int64_t idx[ao_num];
|
||||
double av1[ao_num];
|
||||
for (int64_t k=0 ; k<ao_num ; ++k) {
|
||||
if (avgl1[k] != 0.) {
|
||||
idx[nidx] = k;
|
||||
av1[nidx] = avgl1[k];
|
||||
++nidx;
|
||||
}
|
||||
}
|
||||
|
||||
int64_t n;
|
||||
|
||||
for (n=0 ; n < nidx-4 ; n+=4) {
|
||||
const double* restrict ck1 = coefficient_t + idx[n ]*mo_num;
|
||||
const double* restrict ck2 = coefficient_t + idx[n+1]*mo_num;
|
||||
const double* restrict ck3 = coefficient_t + idx[n+2]*mo_num;
|
||||
const double* restrict ck4 = coefficient_t + idx[n+3]*mo_num;
|
||||
|
||||
const double a11 = av1[n ];
|
||||
const double a21 = av1[n+1];
|
||||
const double a31 = av1[n+2];
|
||||
const double a41 = av1[n+3];
|
||||
|
||||
#ifdef HAVE_OPENMP
|
||||
#pragma omp simd
|
||||
#endif
|
||||
for (int64_t i=0 ; i<mo_num ; ++i) {
|
||||
vgl1[i] = vgl1[i] + ck1[i] * a11 + ck2[i] * a21 + ck3[i] * a31 + ck4[i] * a41;
|
||||
}
|
||||
}
|
||||
|
||||
const int64_t n0 = n < 0 ? 0 : n;
|
||||
for (int64_t m=n0 ; m < nidx ; m+=1) {
|
||||
const double* restrict ck = coefficient_t + idx[m]*mo_num;
|
||||
const double a1 = av1[m];
|
||||
|
||||
#ifdef HAVE_OPENMP
|
||||
#pragma omp simd
|
||||
#endif
|
||||
for (int64_t i=0 ; i<mo_num ; ++i) {
|
||||
vgl1[i] += ck[i] * a1;
|
||||
}
|
||||
}
|
||||
}
|
||||
return QMCKL_SUCCESS;
|
||||
}
|
||||
#endif
|
||||
#+end_src
|
||||
|
||||
** Computation of MOs: values, gradient, Laplacian
|
||||
|
||||
*** Get
|
||||
|
||||
@ -451,7 +920,7 @@ qmckl_get_mo_basis_mo_vgl(qmckl_context context,
|
||||
rc = qmckl_provide_mo_vgl(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
const int64_t sze = ctx->point.num * 5 * ctx->mo_basis.mo_num;
|
||||
@ -507,7 +976,7 @@ qmckl_get_mo_basis_mo_vgl_inplace (qmckl_context context,
|
||||
|
||||
qmckl_exit_code rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
const int64_t sze = ctx->mo_basis.mo_num * 5 * ctx->point.num;
|
||||
@ -563,7 +1032,7 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context)
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!ctx->ao_basis.provided) {
|
||||
@ -639,7 +1108,7 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context)
|
||||
| ~ao_num~ | ~int64_t~ | in | Number of AOs |
|
||||
| ~mo_num~ | ~int64_t~ | in | Number of MOs |
|
||||
| ~point_num~ | ~int64_t~ | in | Number of points |
|
||||
| ~coef_normalized_t~ | ~double[mo_num][ao_num]~ | in | Transpose of the AO to MO transformation matrix |
|
||||
| ~coefficient_t~ | ~double[mo_num][ao_num]~ | in | Transpose of the AO to MO transformation matrix |
|
||||
| ~ao_vgl~ | ~double[point_num][5][ao_num]~ | in | Value, gradients and Laplacian of the AOs |
|
||||
| ~mo_vgl~ | ~double[point_num][5][mo_num]~ | out | Value, gradients and Laplacian of the MOs |
|
||||
|
||||
@ -653,7 +1122,7 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context)
|
||||
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
|
||||
integer function qmckl_compute_mo_basis_mo_vgl_doc_f(context, &
|
||||
ao_num, mo_num, point_num, &
|
||||
coef_normalized_t, ao_vgl, mo_vgl) &
|
||||
coefficient_t, ao_vgl, mo_vgl) &
|
||||
result(info)
|
||||
use qmckl
|
||||
implicit none
|
||||
@ -661,7 +1130,7 @@ integer function qmckl_compute_mo_basis_mo_vgl_doc_f(context, &
|
||||
integer*8 , intent(in) :: ao_num, mo_num
|
||||
integer*8 , intent(in) :: point_num
|
||||
double precision , intent(in) :: ao_vgl(ao_num,5,point_num)
|
||||
double precision , intent(in) :: coef_normalized_t(mo_num,ao_num)
|
||||
double precision , intent(in) :: coefficient_t(mo_num,ao_num)
|
||||
double precision , intent(out) :: mo_vgl(mo_num,5,point_num)
|
||||
integer*8 :: i,j,k
|
||||
double precision :: c1, c2, c3, c4, c5
|
||||
@ -676,15 +1145,21 @@ integer function qmckl_compute_mo_basis_mo_vgl_doc_f(context, &
|
||||
c4 = ao_vgl(k,4,j)
|
||||
c5 = ao_vgl(k,5,j)
|
||||
do i=1,mo_num
|
||||
mo_vgl(i,1,j) = mo_vgl(i,1,j) + coef_normalized_t(i,k) * c1
|
||||
mo_vgl(i,2,j) = mo_vgl(i,2,j) + coef_normalized_t(i,k) * c2
|
||||
mo_vgl(i,3,j) = mo_vgl(i,3,j) + coef_normalized_t(i,k) * c3
|
||||
mo_vgl(i,4,j) = mo_vgl(i,4,j) + coef_normalized_t(i,k) * c4
|
||||
mo_vgl(i,5,j) = mo_vgl(i,5,j) + coef_normalized_t(i,k) * c5
|
||||
mo_vgl(i,1,j) = mo_vgl(i,1,j) + coefficient_t(i,k) * c1
|
||||
mo_vgl(i,2,j) = mo_vgl(i,2,j) + coefficient_t(i,k) * c2
|
||||
mo_vgl(i,3,j) = mo_vgl(i,3,j) + coefficient_t(i,k) * c3
|
||||
mo_vgl(i,4,j) = mo_vgl(i,4,j) + coefficient_t(i,k) * c4
|
||||
mo_vgl(i,5,j) = mo_vgl(i,5,j) + coefficient_t(i,k) * c5
|
||||
end do
|
||||
end if
|
||||
end do
|
||||
end do
|
||||
info = QMCKL_SUCCESS
|
||||
|
||||
! info = qmckl_dgemm(context,'N', 'N', mo_num, point_num, ao_num, 1.d0, &
|
||||
! coefficient_t, int(size(coefficient_t,1),8), &
|
||||
! ao_vgl, int(size(ao_vgl,1),8), 0.d0, &
|
||||
! mo_vgl, int(size(mo_vgl,1),8))
|
||||
|
||||
end function qmckl_compute_mo_basis_mo_vgl_doc_f
|
||||
#+end_src
|
||||
@ -698,7 +1173,7 @@ end function qmckl_compute_mo_basis_mo_vgl_doc_f
|
||||
const int64_t ao_num,
|
||||
const int64_t mo_num,
|
||||
const int64_t point_num,
|
||||
const double* coef_normalized_t,
|
||||
const double* coefficient_t,
|
||||
const double* ao_vgl,
|
||||
double* const mo_vgl );
|
||||
#+end_src
|
||||
@ -712,7 +1187,7 @@ end function qmckl_compute_mo_basis_mo_vgl_doc_f
|
||||
const int64_t ao_num,
|
||||
const int64_t mo_num,
|
||||
const int64_t point_num,
|
||||
const double* coef_normalized_t,
|
||||
const double* coefficient_t,
|
||||
const double* ao_vgl,
|
||||
double* const mo_vgl );
|
||||
#+end_src
|
||||
@ -722,7 +1197,7 @@ end function qmckl_compute_mo_basis_mo_vgl_doc_f
|
||||
#+RESULTS:
|
||||
#+begin_src f90 :tangle (eval f) :comments org :exports none
|
||||
integer(c_int32_t) function qmckl_compute_mo_basis_mo_vgl_doc &
|
||||
(context, ao_num, mo_num, point_num, coef_normalized_t, ao_vgl, mo_vgl) &
|
||||
(context, ao_num, mo_num, point_num, coefficient_t, ao_vgl, mo_vgl) &
|
||||
bind(C) result(info)
|
||||
|
||||
use, intrinsic :: iso_c_binding
|
||||
@ -732,13 +1207,13 @@ end function qmckl_compute_mo_basis_mo_vgl_doc_f
|
||||
integer (c_int64_t) , intent(in) , value :: ao_num
|
||||
integer (c_int64_t) , intent(in) , value :: mo_num
|
||||
integer (c_int64_t) , intent(in) , value :: point_num
|
||||
real (c_double ) , intent(in) :: coef_normalized_t(ao_num,mo_num)
|
||||
real (c_double ) , intent(in) :: coefficient_t(ao_num,mo_num)
|
||||
real (c_double ) , intent(in) :: ao_vgl(ao_num,5,point_num)
|
||||
real (c_double ) , intent(out) :: mo_vgl(mo_num,5,point_num)
|
||||
|
||||
integer(c_int32_t), external :: qmckl_compute_mo_basis_mo_vgl_doc_f
|
||||
info = qmckl_compute_mo_basis_mo_vgl_doc_f &
|
||||
(context, ao_num, mo_num, point_num, coef_normalized_t, ao_vgl, mo_vgl)
|
||||
(context, ao_num, mo_num, point_num, coefficient_t, ao_vgl, mo_vgl)
|
||||
|
||||
end function qmckl_compute_mo_basis_mo_vgl_doc
|
||||
#+end_src
|
||||
@ -749,19 +1224,18 @@ qmckl_compute_mo_basis_mo_vgl (const qmckl_context context,
|
||||
const int64_t ao_num,
|
||||
const int64_t mo_num,
|
||||
const int64_t point_num,
|
||||
const double* coef_normalized_t,
|
||||
const double* coefficient_t,
|
||||
const double* ao_vgl,
|
||||
double* const mo_vgl )
|
||||
{
|
||||
#ifdef HAVE_HPC
|
||||
return qmckl_compute_mo_basis_mo_vgl_hpc (context, ao_num, mo_num, point_num, coef_normalized_t, ao_vgl, mo_vgl);
|
||||
return qmckl_compute_mo_basis_mo_vgl_hpc (context, ao_num, mo_num, point_num, coefficient_t, ao_vgl, mo_vgl);
|
||||
#else
|
||||
return qmckl_compute_mo_basis_mo_vgl_doc (context, ao_num, mo_num, point_num, coef_normalized_t, ao_vgl, mo_vgl);
|
||||
return qmckl_compute_mo_basis_mo_vgl_doc (context, ao_num, mo_num, point_num, coefficient_t, ao_vgl, mo_vgl);
|
||||
#endif
|
||||
}
|
||||
#+end_src
|
||||
|
||||
|
||||
*** HPC version
|
||||
|
||||
|
||||
@ -772,7 +1246,7 @@ qmckl_compute_mo_basis_mo_vgl_hpc (const qmckl_context context,
|
||||
const int64_t ao_num,
|
||||
const int64_t mo_num,
|
||||
const int64_t point_num,
|
||||
const double* coef_normalized_t,
|
||||
const double* coefficient_t,
|
||||
const double* ao_vgl,
|
||||
double* const mo_vgl );
|
||||
#endif
|
||||
@ -785,20 +1259,23 @@ qmckl_compute_mo_basis_mo_vgl_hpc (const qmckl_context context,
|
||||
const int64_t ao_num,
|
||||
const int64_t mo_num,
|
||||
const int64_t point_num,
|
||||
const double* restrict coef_normalized_t,
|
||||
const double* restrict coefficient_t,
|
||||
const double* restrict ao_vgl,
|
||||
double* restrict const mo_vgl )
|
||||
{
|
||||
assert (context != QMCKL_NULL_CONTEXT);
|
||||
|
||||
#ifdef HAVE_OPENMP
|
||||
#pragma omp parallel for
|
||||
#endif
|
||||
for (int64_t ipoint=0 ; ipoint < point_num ; ++ipoint) {
|
||||
double* restrict const vgl1 = &(mo_vgl[ipoint*5*mo_num]);
|
||||
const double* restrict avgl1 = &(ao_vgl[ipoint*5*ao_num]);
|
||||
double* restrict const vgl2 = vgl1 + mo_num;
|
||||
double* restrict const vgl3 = vgl1 + (mo_num << 1);
|
||||
double* restrict const vgl4 = vgl1 + (mo_num << 1) + mo_num;
|
||||
double* restrict const vgl5 = vgl1 + (mo_num << 2);
|
||||
|
||||
const double* restrict avgl1 = &(ao_vgl[ipoint*5*ao_num]);
|
||||
const double* restrict avgl2 = avgl1 + ao_num;
|
||||
const double* restrict avgl3 = avgl1 + (ao_num << 1);
|
||||
const double* restrict avgl4 = avgl1 + (ao_num << 1) + ao_num;
|
||||
@ -820,7 +1297,6 @@ qmckl_compute_mo_basis_mo_vgl_hpc (const qmckl_context context,
|
||||
double av4[ao_num];
|
||||
double av5[ao_num];
|
||||
for (int64_t k=0 ; k<ao_num ; ++k) {
|
||||
const double* restrict ck1 = coef_normalized_t + k*mo_num;
|
||||
if (avgl1[k] != 0.) {
|
||||
idx[nidx] = k;
|
||||
av1[nidx] = avgl1[k];
|
||||
@ -833,12 +1309,12 @@ qmckl_compute_mo_basis_mo_vgl_hpc (const qmckl_context context,
|
||||
}
|
||||
|
||||
int64_t n;
|
||||
|
||||
for (n=0 ; n < nidx-4 ; n+=4) {
|
||||
int64_t k = idx[n];
|
||||
const double* restrict ck1 = coef_normalized_t + idx[n ]*mo_num;
|
||||
const double* restrict ck2 = coef_normalized_t + idx[n+1]*mo_num;
|
||||
const double* restrict ck3 = coef_normalized_t + idx[n+2]*mo_num;
|
||||
const double* restrict ck4 = coef_normalized_t + idx[n+3]*mo_num;
|
||||
const double* restrict ck1 = coefficient_t + idx[n ]*mo_num;
|
||||
const double* restrict ck2 = coefficient_t + idx[n+1]*mo_num;
|
||||
const double* restrict ck3 = coefficient_t + idx[n+2]*mo_num;
|
||||
const double* restrict ck4 = coefficient_t + idx[n+3]*mo_num;
|
||||
|
||||
const double a11 = av1[n ];
|
||||
const double a21 = av1[n+1];
|
||||
@ -877,15 +1353,14 @@ qmckl_compute_mo_basis_mo_vgl_hpc (const qmckl_context context,
|
||||
}
|
||||
}
|
||||
|
||||
int64_t n0 = nidx-4;
|
||||
n0 = n0 < 0 ? 0 : n0;
|
||||
for (int64_t n=n0 ; n < nidx ; n+=1) {
|
||||
const double* restrict ck = coef_normalized_t + idx[n]*mo_num;
|
||||
const double a1 = av1[n];
|
||||
const double a2 = av2[n];
|
||||
const double a3 = av3[n];
|
||||
const double a4 = av4[n];
|
||||
const double a5 = av5[n];
|
||||
const int64_t n0 = n < 0 ? 0 : n;
|
||||
for (int64_t m=n0 ; m < nidx ; m+=1) {
|
||||
const double* restrict ck = coefficient_t + idx[m]*mo_num;
|
||||
const double a1 = av1[m];
|
||||
const double a2 = av2[m];
|
||||
const double a3 = av3[m];
|
||||
const double a4 = av4[m];
|
||||
const double a5 = av5[m];
|
||||
|
||||
#ifdef HAVE_OPENMP
|
||||
#pragma omp simd
|
||||
@ -904,7 +1379,7 @@ qmckl_compute_mo_basis_mo_vgl_hpc (const qmckl_context context,
|
||||
#endif
|
||||
#+end_src
|
||||
|
||||
*** Test
|
||||
** Test
|
||||
|
||||
#+begin_src python :results output :exports none
|
||||
import numpy as np
|
||||
@ -1093,10 +1568,36 @@ assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
assert(qmckl_mo_basis_provided(context));
|
||||
|
||||
rc = qmckl_context_touch(context);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
double mo_value[walk_num*elec_num][chbrclf_mo_num];
|
||||
rc = qmckl_get_mo_basis_mo_value(context, &(mo_value[0][0]), walk_num*elec_num*chbrclf_mo_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
double mo_vgl[walk_num*elec_num][5][chbrclf_mo_num];
|
||||
rc = qmckl_get_mo_basis_mo_vgl(context, &(mo_vgl[0][0][0]), walk_num*elec_num*5*chbrclf_mo_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
for (int i=0 ; i< walk_num*elec_num; ++i) {
|
||||
for (int k=0 ; k< chbrclf_mo_num ; ++k) {
|
||||
assert(fabs(mo_vgl[i][0][k] - mo_value[i][k]) < 1.e-12) ;
|
||||
}
|
||||
}
|
||||
|
||||
rc = qmckl_context_touch(context);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_get_mo_basis_mo_value(context, &(mo_value[0][0]), walk_num*elec_num*chbrclf_mo_num);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
for (int i=0 ; i< walk_num*elec_num; ++i) {
|
||||
for (int k=0 ; k< chbrclf_mo_num ; ++k) {
|
||||
assert(fabs(mo_vgl[i][0][k] - mo_value[i][k]) < 1.e-12) ;
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
// Test overlap of MO
|
||||
//double point_x[10];
|
||||
//double point_y[10];
|
||||
|
@ -125,7 +125,7 @@ qmckl_exit_code qmckl_init_nucleus(qmckl_context context) {
|
||||
return false;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
ctx->nucleus.uninitialized = (1 << 3) - 1;
|
||||
@ -167,7 +167,7 @@ qmckl_get_nucleus_num (const qmckl_context context, int64_t* const num) {
|
||||
"num is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 0;
|
||||
@ -226,7 +226,7 @@ qmckl_get_nucleus_charge (const qmckl_context context,
|
||||
"charge is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 1;
|
||||
@ -293,7 +293,7 @@ qmckl_get_nucleus_rescale_factor (const qmckl_context context,
|
||||
"rescale_factor_kappa is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
assert (ctx->nucleus.rescale_factor_kappa > 0.0);
|
||||
@ -351,7 +351,7 @@ qmckl_get_nucleus_coord (const qmckl_context context,
|
||||
"coord is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int32_t mask = 1 << 2;
|
||||
@ -410,7 +410,7 @@ bool qmckl_nucleus_provided(const qmckl_context context) {
|
||||
return false;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
return ctx->nucleus.provided;
|
||||
@ -422,10 +422,20 @@ bool qmckl_nucleus_provided(const qmckl_context context) {
|
||||
#+NAME:pre2
|
||||
#+begin_src c :exports none
|
||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_NULL_CONTEXT,
|
||||
"qmckl_set_nucleus_*",
|
||||
NULL);
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
|
||||
if (mask != 0 && !(ctx->nucleus.uninitialized & mask)) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_ALREADY_SET,
|
||||
"qmckl_set_nucleus_*",
|
||||
NULL);
|
||||
}
|
||||
#+end_src
|
||||
|
||||
#+NAME:post2
|
||||
@ -452,6 +462,8 @@ qmckl_exit_code
|
||||
qmckl_set_nucleus_num(qmckl_context context,
|
||||
const int64_t num)
|
||||
{
|
||||
int32_t mask = 1 << 0;
|
||||
|
||||
<<pre2>>
|
||||
|
||||
if (num <= 0) {
|
||||
@ -461,8 +473,6 @@ qmckl_set_nucleus_num(qmckl_context context,
|
||||
"num <= 0");
|
||||
}
|
||||
|
||||
int32_t mask = 1 << 0;
|
||||
|
||||
ctx->nucleus.num = num;
|
||||
|
||||
<<post2>>
|
||||
@ -498,6 +508,8 @@ qmckl_set_nucleus_charge(qmckl_context context,
|
||||
const double* charge,
|
||||
const int64_t size_max)
|
||||
{
|
||||
int32_t mask = 1 << 1;
|
||||
|
||||
<<pre2>>
|
||||
|
||||
if (charge == NULL) {
|
||||
@ -510,8 +522,6 @@ qmckl_set_nucleus_charge(qmckl_context context,
|
||||
int64_t num;
|
||||
qmckl_exit_code rc;
|
||||
|
||||
int32_t mask = 1 << 1;
|
||||
|
||||
rc = qmckl_get_nucleus_num(context, &num);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
@ -569,12 +579,12 @@ qmckl_set_nucleus_coord(qmckl_context context,
|
||||
const double* coord,
|
||||
const int64_t size_max)
|
||||
{
|
||||
int32_t mask = 1 << 2;
|
||||
|
||||
<<pre2>>
|
||||
|
||||
qmckl_exit_code rc;
|
||||
|
||||
int32_t mask = 1 << 2;
|
||||
|
||||
const int64_t nucl_num = (int64_t) ctx->nucleus.num;
|
||||
|
||||
if (ctx->nucleus.coord.data != NULL) {
|
||||
@ -641,6 +651,8 @@ qmckl_exit_code
|
||||
qmckl_set_nucleus_rescale_factor(qmckl_context context,
|
||||
const double rescale_factor_kappa)
|
||||
{
|
||||
int32_t mask = 0; // Can be updated
|
||||
|
||||
<<pre2>>
|
||||
|
||||
if (rescale_factor_kappa <= 0.0) {
|
||||
@ -672,7 +684,6 @@ end interface
|
||||
** Test
|
||||
|
||||
#+begin_src c :tangle (eval c_test)
|
||||
const int64_t nucl_num = chbrclf_nucl_num;
|
||||
const double* nucl_charge = chbrclf_charge;
|
||||
const double* nucl_coord = &(chbrclf_nucl_coord[0][0]);
|
||||
const double nucl_rescale_factor_kappa = 2.0;
|
||||
@ -688,13 +699,13 @@ rc = qmckl_get_nucleus_num (context, &n);
|
||||
assert(rc == QMCKL_NOT_PROVIDED);
|
||||
|
||||
|
||||
rc = qmckl_set_nucleus_num (context, nucl_num);
|
||||
rc = qmckl_set_nucleus_num (context, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(!qmckl_nucleus_provided(context));
|
||||
|
||||
rc = qmckl_get_nucleus_num (context, &n);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(n == nucl_num);
|
||||
assert(n == chbrclf_nucl_num);
|
||||
|
||||
double k;
|
||||
rc = qmckl_get_nucleus_rescale_factor (context, &k);
|
||||
@ -709,41 +720,41 @@ rc = qmckl_get_nucleus_rescale_factor (context, &k);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
assert(k == nucl_rescale_factor_kappa);
|
||||
|
||||
double nucl_coord2[3*nucl_num];
|
||||
double nucl_coord2[3*chbrclf_nucl_num];
|
||||
|
||||
rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2, 3*nucl_num);
|
||||
rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2, 3*chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_NOT_PROVIDED);
|
||||
|
||||
rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0]), 3*nucl_num);
|
||||
rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0]), 3*chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
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*chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
for (size_t k=0 ; k<3 ; ++k) {
|
||||
for (int64_t i=0 ; i<nucl_num ; ++i) {
|
||||
assert( nucl_coord[nucl_num*k+i] == nucl_coord2[3*i+k] );
|
||||
for (int64_t i=0 ; i<chbrclf_nucl_num ; ++i) {
|
||||
assert( nucl_coord[chbrclf_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*chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
for (int64_t i=0 ; i<3*nucl_num ; ++i) {
|
||||
for (int64_t i=0 ; i<3*chbrclf_nucl_num ; ++i) {
|
||||
assert( nucl_coord[i] == nucl_coord2[i] );
|
||||
}
|
||||
|
||||
double nucl_charge2[nucl_num];
|
||||
double nucl_charge2[chbrclf_nucl_num];
|
||||
|
||||
rc = qmckl_get_nucleus_charge(context, nucl_charge2, nucl_num);
|
||||
rc = qmckl_get_nucleus_charge(context, nucl_charge2, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_NOT_PROVIDED);
|
||||
|
||||
rc = qmckl_set_nucleus_charge(context, nucl_charge, nucl_num);
|
||||
rc = qmckl_set_nucleus_charge(context, nucl_charge, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_get_nucleus_charge(context, nucl_charge2, nucl_num);
|
||||
rc = qmckl_get_nucleus_charge(context, nucl_charge2, chbrclf_nucl_num);
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
for (int64_t i=0 ; i<nucl_num ; ++i) {
|
||||
for (int64_t i=0 ; i<chbrclf_nucl_num ; ++i) {
|
||||
assert( nucl_charge[i] == nucl_charge2[i] );
|
||||
}
|
||||
assert(qmckl_nucleus_provided(context));
|
||||
@ -784,7 +795,7 @@ qmckl_get_nucleus_nn_distance(qmckl_context context,
|
||||
qmckl_exit_code rc = qmckl_provide_nn_distance(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
const int64_t sze = ctx->nucleus.num * ctx->nucleus.num;
|
||||
@ -828,7 +839,7 @@ qmckl_exit_code qmckl_provide_nn_distance(qmckl_context context)
|
||||
return (char) 0;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!ctx->nucleus.provided) return QMCKL_NOT_PROVIDED;
|
||||
@ -940,10 +951,10 @@ qmckl_exit_code qmckl_compute_nn_distance (
|
||||
|
||||
assert(qmckl_nucleus_provided(context));
|
||||
|
||||
double distance[nucl_num*nucl_num];
|
||||
rc = qmckl_get_nucleus_nn_distance(context, distance, nucl_num*nucl_num);
|
||||
double distance[chbrclf_nucl_num*chbrclf_nucl_num];
|
||||
rc = qmckl_get_nucleus_nn_distance(context, distance, chbrclf_nucl_num*chbrclf_nucl_num);
|
||||
assert(distance[0] == 0.);
|
||||
assert(distance[1] == distance[nucl_num]);
|
||||
assert(distance[1] == distance[chbrclf_nucl_num]);
|
||||
assert(fabs(distance[1]-2.070304721365169) < 1.e-12);
|
||||
|
||||
#+end_src
|
||||
@ -973,7 +984,7 @@ qmckl_get_nucleus_nn_distance_rescaled(qmckl_context context,
|
||||
qmckl_exit_code rc = qmckl_provide_nn_distance_rescaled(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
const int64_t sze = ctx->nucleus.num * ctx->nucleus.num;
|
||||
@ -1019,7 +1030,7 @@ qmckl_exit_code qmckl_provide_nn_distance_rescaled(qmckl_context context)
|
||||
return (char) 0;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!ctx->nucleus.provided) return QMCKL_NOT_PROVIDED;
|
||||
@ -1167,7 +1178,7 @@ qmckl_exit_code qmckl_get_nucleus_repulsion(qmckl_context context, double* const
|
||||
qmckl_exit_code rc = qmckl_provide_nucleus_repulsion(context);
|
||||
if (rc != QMCKL_SUCCESS) return rc;
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
*energy = ctx->nucleus.repulsion;
|
||||
@ -1203,7 +1214,7 @@ qmckl_exit_code qmckl_provide_nucleus_repulsion(qmckl_context context)
|
||||
return (char) 0;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
qmckl_exit_code rc;
|
||||
|
@ -141,7 +141,7 @@ qmckl_exit_code qmckl_set_numprec_precision(const qmckl_context context, const i
|
||||
"precision > 53");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
|
||||
/* This should be always true because the context is valid */
|
||||
assert (ctx != NULL);
|
||||
@ -185,7 +185,7 @@ int qmckl_get_numprec_precision(const qmckl_context context) {
|
||||
"");
|
||||
}
|
||||
|
||||
const qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
const qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
return ctx->numprec.precision;
|
||||
}
|
||||
#+end_src
|
||||
@ -232,7 +232,7 @@ qmckl_exit_code qmckl_set_numprec_range(const qmckl_context context, const int r
|
||||
"range > 11");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
|
||||
/* This should be always true because the context is valid */
|
||||
assert (ctx != NULL);
|
||||
@ -275,7 +275,7 @@ int qmckl_get_numprec_range(const qmckl_context context) {
|
||||
"");
|
||||
}
|
||||
|
||||
const qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
const qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
return ctx->numprec.range;
|
||||
}
|
||||
#+end_src
|
||||
|
@ -78,10 +78,12 @@ int main() {
|
||||
The following data stored in the context:
|
||||
|
||||
| Variable | Type | Description |
|
||||
|----------+----------------+-------------------------------------------|
|
||||
|--------------+----------------+-------------------------------------------|
|
||||
| ~num~ | ~int64_t~ | Total number of points |
|
||||
| ~alloc_num~ | ~int64_t~ | Numer of allocated number of points |
|
||||
| ~date~ | ~uint64_t~ | Last modification date of the coordinates |
|
||||
| ~coord~ | ~qmckl_matrix~ | ~num~ \times 3 matrix3 |
|
||||
| ~alloc_date~ | ~uint64_t~ | Last modification date of the allocation |
|
||||
| ~coord~ | ~qmckl_matrix~ | ~num~ \times 3 matrix |
|
||||
|
||||
We consider that the matrix is stored 'transposed' and 'normal'
|
||||
corresponds to the 3 \times ~num~ matrix.
|
||||
@ -91,7 +93,9 @@ int main() {
|
||||
#+begin_src c :comments org :tangle (eval h_private_type)
|
||||
typedef struct qmckl_point_struct {
|
||||
int64_t num;
|
||||
int64_t alloc_num;
|
||||
uint64_t date;
|
||||
uint64_t alloc_date;
|
||||
qmckl_matrix coord;
|
||||
} qmckl_point_struct;
|
||||
|
||||
@ -108,7 +112,7 @@ qmckl_exit_code qmckl_init_point(qmckl_context context) {
|
||||
return false;
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
memset(&(ctx->point), 0, sizeof(qmckl_point_struct));
|
||||
@ -148,7 +152,7 @@ qmckl_get_point_num (const qmckl_context context, int64_t* const num) {
|
||||
"num is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
assert (ctx->point.num >= (int64_t) 0);
|
||||
@ -202,7 +206,7 @@ qmckl_get_point(const qmckl_context context,
|
||||
"coord is a null pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
int64_t point_num = ctx->point.num;
|
||||
@ -263,8 +267,9 @@ end interface
|
||||
#+begin_src c :comments org :tangle (eval h_func)
|
||||
qmckl_exit_code qmckl_set_point (qmckl_context context,
|
||||
const char transp,
|
||||
const int64_t num,
|
||||
const double* coord,
|
||||
const int64_t num);
|
||||
const int64_t size_max);
|
||||
#+end_src
|
||||
|
||||
Copy a sequence of ~num~ points $(x,y,z)$ into the context.
|
||||
@ -273,14 +278,22 @@ qmckl_exit_code qmckl_set_point (qmckl_context context,
|
||||
qmckl_exit_code
|
||||
qmckl_set_point (qmckl_context context,
|
||||
const char transp,
|
||||
const int64_t num,
|
||||
const double* coord,
|
||||
const int64_t num)
|
||||
const int64_t size_max)
|
||||
{
|
||||
|
||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||
return QMCKL_NULL_CONTEXT;
|
||||
}
|
||||
|
||||
if (size_max < 3*num) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_INVALID_ARG_4,
|
||||
"qmckl_set_point",
|
||||
"Array too small");
|
||||
}
|
||||
|
||||
if (transp != 'N' && transp != 'T') {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_INVALID_ARG_2,
|
||||
@ -295,11 +308,11 @@ qmckl_set_point (qmckl_context context,
|
||||
"coord is a NULL pointer");
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
qmckl_exit_code rc;
|
||||
if (ctx->point.num < num) {
|
||||
if (num > ctx->point.alloc_num) {
|
||||
|
||||
if (ctx->point.coord.data != NULL) {
|
||||
rc = qmckl_matrix_free(context, &(ctx->point.coord));
|
||||
@ -313,7 +326,6 @@ qmckl_set_point (qmckl_context context,
|
||||
"qmckl_set_point",
|
||||
NULL);
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
ctx->point.num = num;
|
||||
@ -341,6 +353,11 @@ qmckl_set_point (qmckl_context context,
|
||||
rc = qmckl_context_touch(context);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
|
||||
if (num > ctx->point.alloc_num) {
|
||||
ctx->point.alloc_num = num;
|
||||
ctx->point.alloc_date = ctx->point.date;
|
||||
};
|
||||
|
||||
return QMCKL_SUCCESS;
|
||||
|
||||
}
|
||||
@ -349,15 +366,16 @@ qmckl_set_point (qmckl_context context,
|
||||
#+begin_src f90 :comments org :tangle (eval fh_func) :noweb yes
|
||||
interface
|
||||
integer(c_int32_t) function qmckl_set_point(context, &
|
||||
transp, coord, num) bind(C)
|
||||
transp, num, coord, size_max) bind(C)
|
||||
use, intrinsic :: iso_c_binding
|
||||
import
|
||||
implicit none
|
||||
|
||||
integer (c_int64_t) , intent(in) , value :: context
|
||||
character(c_char) , intent(in) , value :: transp
|
||||
real (c_double ) , intent(in) :: coord(*)
|
||||
integer (c_int64_t) , intent(in) , value :: num
|
||||
real (c_double ) , intent(in) :: coord(*)
|
||||
integer (c_int64_t) , intent(in) , value :: size_max
|
||||
end function
|
||||
end interface
|
||||
#+end_src
|
||||
@ -380,7 +398,7 @@ double coord3[point_num*3];
|
||||
rc = qmckl_get_point (context, 'N', coord2, (point_num*3));
|
||||
assert(rc == QMCKL_NOT_PROVIDED);
|
||||
|
||||
rc = qmckl_set_point (context, 'N', coord, point_num);
|
||||
rc = qmckl_set_point (context, 'N', point_num, coord, (point_num*3));
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
int64_t n;
|
||||
@ -404,7 +422,7 @@ for (int64_t i=0 ; i<point_num ; ++i) {
|
||||
assert( coord[3*i+2] == coord2[i+point_num*2] );
|
||||
}
|
||||
|
||||
rc = qmckl_set_point (context, 'T', coord2, point_num);
|
||||
rc = qmckl_set_point (context, 'T', point_num, coord2, (point_num*3));
|
||||
assert(rc == QMCKL_SUCCESS);
|
||||
|
||||
rc = qmckl_get_point (context, 'N', coord3, (point_num*3));
|
||||
|
@ -965,7 +965,7 @@ qmckl_exit_code qmckl_sherman_morrison_smw32s(const qmckl_context context,
|
||||
rc = qmckl_woodbury_3(context, LDS, Dim, Updates_3block, Updates_index_3block, breakdown, Slater_inv, determinant);
|
||||
if (rc != 0) { // Send the entire block to slagel_splitting
|
||||
uint64_t l = 0;
|
||||
rc = qmckl_slagel_splitting(LDS, Dim, 3, Updates_3block, Updates_index_3block,
|
||||
(void) qmckl_slagel_splitting(LDS, Dim, 3, Updates_3block, Updates_index_3block,
|
||||
breakdown, Slater_inv, later_updates + (Dim * later), later_index + later, &l, determinant);
|
||||
later = later + l;
|
||||
}
|
||||
|
@ -428,13 +428,6 @@ qmckl_trexio_read_ao_X(qmckl_context context, trexio_t* const file)
|
||||
}
|
||||
|
||||
/* Reformat data */
|
||||
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) {
|
||||
const int k = tmp_array[i];
|
||||
if (k < 0 || k >= nucleus_num) {
|
||||
@ -1086,7 +1079,7 @@ qmckl_trexio_read(const qmckl_context context, const char* file_name, const int6
|
||||
|
||||
qmckl_exit_code rc;
|
||||
char file_name_new[size_max+1];
|
||||
strncpy(file_name_new, file_name, size_max+1);
|
||||
strncpy(file_name_new, file_name, size_max);
|
||||
file_name_new[size_max] = '\0';
|
||||
|
||||
#ifdef HAVE_TREXIO
|
||||
|
@ -18,3 +18,4 @@ qmckl_utils.org
|
||||
qmckl_trexio.org
|
||||
qmckl_tests.org
|
||||
qmckl_verificarlo.org
|
||||
qmckl_examples.org
|
||||
|
@ -49,6 +49,10 @@ import_array();
|
||||
/* Include typemaps generated by the process_header.py script */
|
||||
%include "qmckl_include.i"
|
||||
|
||||
/* Handle properly get_point */
|
||||
|
||||
|
||||
|
||||
/* exception.i is a generic (language-independent) module */
|
||||
%include "exception.i"
|
||||
|
||||
|
159
tools/lib.org
159
tools/lib.org
@ -15,7 +15,6 @@
|
||||
(org-entry-get nil key t))
|
||||
#+end_src
|
||||
|
||||
|
||||
** Table of function arguments
|
||||
|
||||
#+NAME: test
|
||||
@ -149,7 +148,6 @@ return template
|
||||
const int64_t ldc );
|
||||
#+end_src
|
||||
|
||||
|
||||
*** Generates a C interface to the Fortran function
|
||||
|
||||
#+NAME: generate_c_interface
|
||||
@ -259,3 +257,160 @@ return results
|
||||
|
||||
|
||||
|
||||
** Creating provide functions
|
||||
|
||||
#+NAME: write_provider_header
|
||||
#+BEGIN_SRC python :var group="GROUP" :var data="DATA" :results drawer :noweb yes :wrap "src c :comments org :tangle (eval h_private_func) :noweb yes :export none"
|
||||
template = "qmckl_exit_code qmckl_provide_{{ group }}_{{ data }}(qmckl_context context);"
|
||||
|
||||
msg = template.replace("{{ group }}", group) \
|
||||
.replace("{{ data }}", data)
|
||||
return msg
|
||||
#+END_SRC
|
||||
|
||||
#+RESULTS: write_provider_header
|
||||
#+begin_src c :comments org :tangle (eval h_private_func) :noweb yes :export none
|
||||
qmckl_exit_code qmckl_provide_GROUP_DATA(qmckl_context context);
|
||||
#+end_src
|
||||
|
||||
#+NAME: write_provider_pre
|
||||
#+BEGIN_SRC python :var group="GROUP" :var data="DATA" :var dimension="DIMENSION" :results drawer :noweb yes :wrap "src c :comments org :tangle (eval c) :noweb yes :export none"
|
||||
template = """qmckl_exit_code qmckl_provide_{{ group }}_{{ data }}(qmckl_context context)
|
||||
{
|
||||
|
||||
qmckl_exit_code rc;
|
||||
|
||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_INVALID_CONTEXT,
|
||||
"qmckl_provide_{{ group }}_{{ data }}",
|
||||
NULL);
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!ctx->{{ group }}.provided) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_NOT_PROVIDED,
|
||||
"qmckl_provide_{{ group }}_{{ data }}",
|
||||
NULL);
|
||||
}
|
||||
|
||||
/* Compute if necessary */
|
||||
if (ctx->point.date > ctx->{{ group }}.{{ data }}_date) {
|
||||
|
||||
if (ctx->point.alloc_date > ctx->{{ group }}.{{ data }}_date) {
|
||||
if (ctx->{{ group }}.{{ data }} != NULL) {
|
||||
rc = qmckl_free(context, ctx->{{ group }}.{{ data }});
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
ctx->{{ group }}.{{ data }} = NULL;
|
||||
}
|
||||
}
|
||||
|
||||
/* Allocate array */
|
||||
if (ctx->{{ group }}.{{ data }} == NULL) {
|
||||
|
||||
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
|
||||
mem_info.size = {{ dimension }} * sizeof(double);
|
||||
double* {{ data }} = (double*) qmckl_malloc(context, mem_info);
|
||||
|
||||
if ({{ data }} == NULL) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_ALLOCATION_FAILED,
|
||||
"qmckl_{{ group }}_{{ data }}",
|
||||
NULL);
|
||||
}
|
||||
ctx->{{ group }}.{{ data }} = {{ data }};
|
||||
}
|
||||
|
||||
"""
|
||||
|
||||
msg = template.replace("{{ group }}", group) \
|
||||
.replace("{{ data }}", data) \
|
||||
.replace("{{ dimension }}", dimension)
|
||||
return msg
|
||||
#+END_SRC
|
||||
|
||||
#+RESULTS: write_provider_pre
|
||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :export none
|
||||
qmckl_exit_code qmckl_provide_GROUP_DATA(qmckl_context context)
|
||||
{
|
||||
|
||||
qmckl_exit_code rc;
|
||||
|
||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_INVALID_CONTEXT,
|
||||
"qmckl_provide_GROUP_DATA",
|
||||
NULL);
|
||||
}
|
||||
|
||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||
assert (ctx != NULL);
|
||||
|
||||
if (!ctx->GROUP.provided) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_NOT_PROVIDED,
|
||||
"qmckl_provide_GROUP_DATA",
|
||||
NULL);
|
||||
}
|
||||
|
||||
/* Compute if necessary */
|
||||
if (ctx->point.date > ctx->GROUP.DATA_date) {
|
||||
|
||||
if (ctx->point.alloc_date > ctx->GROUP.DATA_date) {
|
||||
rc = qmckl_free(context, ctx->GROUP.DATA);
|
||||
assert (rc == QMCKL_SUCCESS);
|
||||
ctx->GROUP.DATA = NULL;
|
||||
}
|
||||
|
||||
/* Allocate array */
|
||||
if (ctx->GROUP.DATA == NULL) {
|
||||
|
||||
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
|
||||
mem_info.size = DIMENSION * sizeof(double);
|
||||
double* DATA = (double*) qmckl_malloc(context, mem_info);
|
||||
|
||||
if (DATA == NULL) {
|
||||
return qmckl_failwith( context,
|
||||
QMCKL_ALLOCATION_FAILED,
|
||||
"qmckl_GROUP_DATA",
|
||||
NULL);
|
||||
}
|
||||
ctx->GROUP.DATA = DATA;
|
||||
}
|
||||
|
||||
#+end_src
|
||||
|
||||
#+NAME: write_provider_post
|
||||
#+BEGIN_SRC python :var group="BASIS" :var data="DATA" :results drawer :noweb yes :wrap "src c :comments org :tangle (eval c) :noweb yes :export none"
|
||||
template = """ if (rc != QMCKL_SUCCESS) {
|
||||
return rc;
|
||||
}
|
||||
|
||||
ctx->{{ group }}.{{ data }}_date = ctx->date;
|
||||
}
|
||||
|
||||
return QMCKL_SUCCESS;
|
||||
}
|
||||
"""
|
||||
|
||||
msg = template.replace("{{ group }}", group) \
|
||||
.replace("{{ data }}", data)
|
||||
|
||||
return msg
|
||||
#+END_SRC
|
||||
|
||||
#+RESULTS: write_provider_post
|
||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :export none
|
||||
if (rc != QMCKL_SUCCESS) {
|
||||
return rc;
|
||||
}
|
||||
|
||||
ctx->BASIS.DATA_date = ctx->date;
|
||||
}
|
||||
|
||||
return QMCKL_SUCCESS;
|
||||
}
|
||||
#+end_src
|
||||
|
Loading…
Reference in New Issue
Block a user