From abfce6b4be92a35e77d7fe856c7345167548d06f Mon Sep 17 00:00:00 2001 From: scemama Date: Fri, 19 Mar 2021 12:49:11 +0000 Subject: [PATCH] =?UTF-8?q?Deploying=20to=20gh-pages=20from=20@=20TREX-CoE?= =?UTF-8?q?/qmckl@87ad36e342996f0a19429e5cb78933b7d5e472ae=20=F0=9F=9A=80?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- README.html | 355 +++++ config.el | 98 -- index.html | 3090 ++----------------------------------------- qmckl.css | 972 ++++++++++++++ qmckl.html | 720 ++++++++++ qmckl_ao.html | 1206 +++++++++++++++++ qmckl_context.html | 1869 ++++++++++++++++++++++++++ qmckl_distance.html | 628 +++++++++ qmckl_error.html | 452 +++++++ qmckl_memory.html | 436 ++++++ test_qmckl.html | 318 +++++ theme.setup | 15 + 12 files changed, 7078 insertions(+), 3081 deletions(-) create mode 100644 README.html delete mode 100755 config.el create mode 100644 qmckl.css create mode 100644 qmckl.html create mode 100644 qmckl_ao.html create mode 100644 qmckl_context.html create mode 100644 qmckl_distance.html create mode 100644 qmckl_error.html create mode 100644 qmckl_memory.html create mode 100644 test_qmckl.html create mode 100644 theme.setup diff --git a/README.html b/README.html new file mode 100644 index 0000000..27c961a --- /dev/null +++ b/README.html @@ -0,0 +1,355 @@ + + + + + + + +QMCkl source code documentation + + + + + + + + + + + +
+ UP + | + HOME +
+

QMCkl source code documentation

+
+ + + + + +
+ + +

+The ultimate goal of the QMCkl library is to provide a high-performance +implementation of the main kernels of QMC. In this particular +implementation of the library, we focus on the definition of the API +and the tests, and on a pedagogical presentation of the +algorithms. We expect the HPC experts to use this repository as a +reference for re-writing optimized libraries. +

+ +

+The source code of the library is available at +https://github.com/trex-coe/qmckl +and bug reports should be submitted at +https://github.com/trex-coe/qmckl/issues. +

+ +
+ +

+euflag.jpg TREX: Targeting Real Chemical Accuracy at the Exascale project has received funding from the European Union’s Horizon 2020 - Research and Innovation program - under grant agreement no. 952165. The content of this document does not represent the opinion of the European Union, and the European Union is not responsible for any use that might be made of such content. +

+
+
+

Author: TREX CoE

+

Created: 2021-03-19 Fri 12:49

+

Validate

+
+ + diff --git a/config.el b/config.el deleted file mode 100755 index 0a7f879..0000000 --- a/config.el +++ /dev/null @@ -1,98 +0,0 @@ -;; Thanks to Tobias's answer on Emacs Stack Exchange: -;; https://emacs.stackexchange.com/questions/38437/org-mode-batch-export-missing-syntax-highlighting - -(package-initialize) -(add-to-list 'package-archives - '("gnu" . "https://elpa.gnu.org/packages/")) -(add-to-list 'package-archives - '("melpa-stable" . "https://stable.melpa.org/packages/")) -(add-to-list 'package-archives - '("melpa" . "https://melpa.org/packages/")) -(setq package-archive-priorities '(("melpa-stable" . 100) - ("melpa" . 50) - ("gnu" . 10))) - - -(require 'htmlize) -(require 'font-lock) -(require 'subr-x) ;; for `when-let' -(setq org-confirm-babel-evaluate nil) -(global-font-lock-mode t) - -(org-babel-do-load-languages - 'org-babel-load-languages - '( - (emacs-lisp . t) - (shell . t) - (python . t) - (C . t) - (org . t) - (makefile . t) - )) - - - - -(unless (boundp 'maximal-integer) - (defconst maximal-integer (lsh -1 -1) - "Maximal integer value representable natively in emacs lisp.")) - -(defun face-spec-default (spec) - "Get list containing at most the default entry of face SPEC. -Return nil if SPEC has no default entry." - (let* ((first (car-safe spec)) - (display (car-safe first))) - (when (eq display 'default) - (list (car-safe spec))))) - -(defun face-spec-min-color (display-atts) - "Get min-color entry of DISPLAY-ATTS pair from face spec." - (let* ((display (car-safe display-atts))) - (or (car-safe (cdr (assoc 'min-colors display))) - maximal-integer))) - -(defun face-spec-highest-color (spec) - "Search face SPEC for highest color. -That means the DISPLAY entry of SPEC -with class 'color and highest min-color value." - (let ((color-list (cl-remove-if-not - (lambda (display-atts) - (when-let ((display (car-safe display-atts)) - (class (and (listp display) - (assoc 'class display))) - (background (assoc 'background display))) - (and (member 'light (cdr background)) - (member 'color (cdr class))))) - spec))) - (cl-reduce (lambda (display-atts1 display-atts2) - (if (> (face-spec-min-color display-atts1) - (face-spec-min-color display-atts2)) - display-atts1 - display-atts2)) - (cdr color-list) - :initial-value (car color-list)))) - -(defun face-spec-t (spec) - "Search face SPEC for fall back." - (cl-find-if (lambda (display-atts) - (eq (car-safe display-atts) t)) - spec)) - -(defun my-face-attribute (face attribute &optional frame inherit) - "Get FACE ATTRIBUTE from `face-user-default-spec' and not from `face-attribute'." - (let* ((face-spec (face-user-default-spec face)) - (display-attr (or (face-spec-highest-color face-spec) - (face-spec-t face-spec))) - (attr (cdr display-attr)) - (val (or (plist-get attr attribute) (car-safe (cdr (assoc attribute attr)))))) - ;; (message "attribute: %S" attribute) ;; for debugging - (when (and (null (eq attribute :inherit)) - (null val)) - (let ((inherited-face (my-face-attribute face :inherit))) - (when (and inherited-face - (null (eq inherited-face 'unspecified))) - (setq val (my-face-attribute inherited-face attribute))))) - ;; (message "face: %S attribute: %S display-attr: %S, val: %S" face attribute display-attr val) ;; for debugging - (or val 'unspecified))) - -(advice-add 'face-attribute :override #'my-face-attribute) diff --git a/index.html b/index.html index 28222ab..27c961a 100644 --- a/index.html +++ b/index.html @@ -3,11 +3,12 @@ "http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd"> - + QMCkl source code documentation + - - - - - - + + + + + - - -
+
+ UP + | + HOME +

QMCkl source code documentation

- +
+ + + + + +
-
-

1 Introduction

-

The ultimate goal of the QMCkl library is to provide a high-performance implementation of the main kernels of QMC. In this particular -implementation of the library, we focus on the definition of the API and the tests, and -on a pedagogical presentation of the algorithms. We expect the -HPC experts to use this repository as a reference for re-writing -optimized libraries. -

-
- - -
-

1.1 Literate programming

-
-

-In a traditional source code, most of the lines of source files of a program -are code, scripts, Makefiles, and only a few lines are comments explaining -parts of the code that are non-trivial to understand. The documentation of -the prorgam is usually written in a separate directory, and is often outdated -compared to the code. -

- -

-Literate programming is a different approach to programming, -where the program is considered as a publishable-quality document. Most of -the lines of the source files are text, mathematical formulas, tables, -figures, etc, and the lines of code are just the translation in a computer -language of the ideas and algorithms expressed in the text. More importantly, -the "document" is structured like a text document with sections, subsections, -a bibliography, a table of contents etc, and the place where pieces of code -appear are the places where they should belong for the reader to understand -the logic of the program, not the places where the compiler expects to find -them. Both the publishable-quality document and the binary executable are -produced from the same source files. -

- -

-Literate programming is particularly well adapted in this context, as the -central part of this project is the documentation of an API. The -implementation of the algorithms is just an expression of the algorithms in a -language that can be compiled, so that the correctness of the algorithms can -be tested. -

- -

-We have chosen to write the source files in org-mode format, -as any text editor can be used to edit org-mode files. To -produce the documentation, there exists multiple possibilities to convert -org-mode files into different formats such as HTML or PDF. The source code is -easily extracted from the org-mode files invoking the Emacs text editor from -the command-line in the Makefile, and then the produced files are compiled. -Moreover, within the Emacs text editor the source code blocks can be executed -interactively, in the same spirit as Jupyter notebooks. -

-
-
- - -
-

1.2 Source code editing

-
-

-For a tutorial on literate programming with org-mode, follow this link. -

- -

-Any text editor can be used to edit org-mode files. For a better -user experience Emacs is recommended. For users hating Emacs, it -is good to know that Emacs can behave like Vim when switched into -``Evil'' mode. -

- -

-In the tools/init.el file, we provide a minimal Emacs configuration -file for vim users. This file should be copied into .emacs.d/init.el. -

- -

-For users with a preference for Jupyter notebooks, we also provide the -tools/nb_to_org.sh script can convert jupyter notebooks into org-mode -files. -

- -

-Note that pandoc can be used to convert multiple markdown formats into -org-mode. -

-
-
- - -
-

1.3 Choice of the programming language

-
-

-Most of the codes of the TREX CoE are written in Fortran with some scripts in -Bash and Python. Outside of the CoE, Fortran is also important (Casino, Amolqc), -and other important languages used by the community are C and C++ (QMCPack, -QWalk), and Julia is gaining in popularity. The library we design should be -compatible with all of these languages. The QMCkl API has to be compatible -with the C language since libraries with a C-compatible API can be used in -every other language. -

- -

-High-performance versions of the QMCkl, with the same API, will be rewritten by -the experts in HPC. These optimized libraries will be tuned for specific -architectures, among which we can cite x86 based processors, and GPU -accelerators. Nowadays, the most efficient software tools to take advantage of -low-level features of the processor (intrinsics) and of GPUs are for C++ -developers. It is highly probable that the optimized implementations will be -written in C++, and this is agreement with our choice to make the API -C-compatible. -

- -

-Fortran is one of the most common languages used by the community, and is simple -enough to make the algorithms readable both by experts in QMC, and experts in -HPC. Hence we propose in this pedagogical implementation of QMCkl to use Fortran -to express the QMC algorithms. As the main languages of the library is C, this -implies that the exposed C functions call the Fortran routine. However, for -internal functions related to system programming, the C language is more natural -than Fortran. -

- -

-The Fortran source files should provide a C interface using the -iso_c_binding module. The name of the Fortran source files should end with -_f.f90 to be properly handled by the Makefile. The names of the functions -defined in Fortran should be the same as those exposed in the API suffixed by -_f. Fortran interfaces should also be written in the qmckl_f.f90 file. -

- -

-For more guidelines on using Fortran to generate a C interface, see -this link. -

-
-
- -
-

1.4 Design of the library

-
-

-The proposed API should allow the library to: deal with memory transfers -between CPU and accelerators, and to use different levels of floating-point -precision. We chose a multi-layered design with low-level and high-level -functions (see below). -

-
- -
-

1.4.1 Naming conventions

-
-

-To avoid namespace collisions, we use qmckl_ as a prefix for all exported -functions and variables. All exported header files should have a file name -prefixed with qmckl_. -

- -

-If the name of the org-mode file is xxx.org, the name of the -produced C files should be xxx.c and xxx.h and the name of the -produced Fortran file should be xxx.f90. -

- -

-Arrays are in uppercase and scalars are in lowercase. -

- -

-In the names of the variables and functions, only the singular -form is allowed. -

-
-
- -
-

1.4.2 Application programming interface

-
-

-In the C language, the number of bits used by the integer types can change -from one architecture to another one. To circumvent this problem, we choose to -use the integer types defined in <stdint.h> where the number of bits used for -the integers are fixed. -

- -

-To ensure that the library will be easily usable in any other language -than C, we restrict the data types in the interfaces to the following: -

-
    -
  • 32-bit and 64-bit integers, scalars and and arrays (int32_t and int64_t)
    • -
  • 32-bit and 64-bit floats, scalars and and arrays (float and double)
    • -
  • Pointers are always casted into 64-bit integers, even on legacy 32-bit architectures
    • -
  • ASCII strings are represented as a pointers to character arrays -and terminated by a '\0' character (C convention).
    • -
  • Complex numbers can be represented by an array of 2 floats.
    • -
  • Boolean variables are stored as integers, 1 for true and 0 for false
    • -
  • Floating point variables should be by default
    • -
  • double unless explicitly mentioned
    • -
  • integers used for counting should always be int64_t
    • -
- -

-To facilitate the use in other languages than C, we will provide some -bindings in other languages in other repositories. -

-
-
- -
-

1.4.3 Global state

-
-

-Global variables should be avoided in the library, because it is -possible that one single program needs to use multiple instances -of the library. To solve this problem we propose to use a pointer -to a context variable, built by the library with the -qmckl_context_create function. The context contains the global -state of the library, and is used as the first argument of many -QMCkl functions. -

- -

-The internal structure of the context is not specified, to give a -maximum of freedom to the different implementations. Modifying -the state is done by setters and getters, prefixed by -qmckl_context_set_ an qmckl_context_get_. When a context -variable is modified by a setter, a copy of the old data structure -is made and updated, and the pointer to the new data structure is -returned, such that the old contexts can still be accessed. It is -also possible to modify the state in an impure fashion, using the -qmckl_context_update_ functions. The context and its old -versions can be destroyed with qmckl_context_destroy. -

-
-
- -
-

1.4.4 Low-level functions

-
-

-Low-level functions are very simple functions which are leaves of -the function call tree (they don't call any other QMCkl function). -

- -

-These functions are pure, and unaware of the QMCkl -context. They are not allowed to allocate/deallocate memory, and -if they need temporary memory it should be provided in input. -

-
-
- -
-

1.4.5 High-level functions

-
-

-High-level functions are at the top of the function call tree. -They are able to choose which lower-level function to call -depending on the required precision, and do the corresponding type -conversions. These functions are also responsible for allocating -temporary storage, to simplify the use of accelerators. -

- -

-The high-level functions should be pure, unless the introduction -of non-purity is justified. All the side effects should be made in -the context variable. -

-
-
- -
-

1.4.6 Numerical precision

-
-

-The number of bits of precision required for a function should be -given as an input of low-level computational functions. This input -will be used to define the values of the different thresholds that -might be used to avoid computing unnecessary noise. High-level -functions will use the precision specified in the context -variable. -

-
-
-
- -
-

1.5 Algorithms

-
-

-Reducing the scaling of an algorithm usually implies also reducing -its arithmetic complexity (number of flops per byte). Therefore, -for small sizes \(\mathcal{O}(N^3)\) and \(\mathcal{O}(N^2)\) -algorithms are better adapted than linear scaling algorithms. As -QMCkl is a general purpose library, multiple algorithms should be -implemented adapted to different problem sizes. -

-
-
- -
-

1.6 Rules for the API

-
-
    -
  • stdint should be used for integers (int32_t, int64_t)
    • -
  • integers used for counting should always be int64_t
    • -
  • floats should be by default double, unless explicitly mentioned
    • -
  • pointers are converted to int64_t to increase portability
    • -
-
-
-
- -
-

2 Documentation

-
-
- - -
-

2.1 qmckl.h header file

-
-

-The qmckl.h header file has to be included in C codes when -QMCkl functions are used: -#+BEGINSRC C :tangle none -#include "qmckl.h" -#+ENDSRC f90 -

- - -

-In Fortran programs, the qmckl_f.f90 interface file should be -included in the source code using the library, and the Fortran codes -should use the qmckl module as -#+BEGINSRC f90 :tangle none -use qmckl -#+ENDSRC f90 -

-
-
- -
-

2.2 Context

-
-

-This file is written in C because it is more natural to express the -context in C than in Fortran. -

- -

-2 files are produced: -

-
    -
  • a source file : qmckl_context.c
    • -
  • a test file : test_qmckl_context.c
    • -
-
- -
-

2.2.1 Context

-
-

-The context variable is a handle for the state of the library, and -is stored in the following data structure, which can't be seen -outside of the library. To simplify compatibility with other -languages, the pointer to the internal data structure is converted -into a 64-bit signed integer, defined in the qmckl_context type. -A value of 0 for the context is equivalent to a NULL pointer. -

- -
-
typedef int64_t qmckl_context ;
-
-
-
- -
    -
  1. Data for error handling
    -
    -

    -We define here the the data structure containing the strings -necessary for error handling. -

    - -
    -
    #define  QMCKL_MAX_FUN_LEN  256
-#define  QMCKL_MAX_MSG_LEN  1024
-
-typedef struct qmckl_error_struct {
-
-  qmckl_exit_code exit_code;
-  char function[QMCKL_MAX_FUN_LEN];
-  char message [QMCKL_MAX_MSG_LEN];
-
-} qmckl_error_struct;
-
    -
    -
    -
    2. - - -
  2. Basis set data structure
    -
    -

    -Data structure for the info related to the atomic orbitals -basis set. -

    - -
    -
    typedef struct qmckl_ao_basis_struct {
-
-  int64_t   shell_num;
-  int64_t   prim_num;
-  int64_t * shell_center;
-  int32_t * shell_ang_mom;
-  double  * shell_factor;
-  double  * exponent    ;
-  double  * coefficient ;
-  int64_t * shell_prim_num;
-  char      type;
-
-} qmckl_ao_basis_struct;
-
    -
    -
    - -
      -
    1. Source
      -
      -

      -The tag is used internally to check if the memory domain pointed -by a pointer is a valid context. -

      - -
      -
      typedef struct qmckl_context_struct {
-
-  struct qmckl_context_struct * prev;
-
-  /* Molecular system */
-  //  struct qmckl_nucleus_struct     * nucleus;
-  //  struct qmckl_electron_struct    * electron;
-  struct qmckl_ao_basis_struct      * ao_basis;
-  //  struct qmckl_mo_struct          * mo;
-  //  struct qmckl_determinant_struct * det;
-
-  /* Numerical precision */
-  uint32_t tag;
-  int32_t precision;
-  int32_t range;
-
-  /* Error handling */
-  struct qmckl_error_struct   * error;
-
-} qmckl_context_struct;
-
-#define VALID_TAG   0xBEEFFACE
-#define INVALID_TAG 0xDEADBEEF
-
      -
      -
      -
      2. -
    -
    3. - -
  3. qmckl_context_update_error
    -
    -
    -
    qmckl_exit_code
-qmckl_context_update_error(qmckl_context context, const qmckl_exit_code exit_code, const char* function, const char* message);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_exit_code
-qmckl_context_update_error(qmckl_context context, const qmckl_exit_code exit_code, const char* function, const char* message)
-{
-  assert (context != 0);
-  assert (function != NULL);
-  assert (message != NULL);
-  assert (exit_code > 0);
-  assert (exit_code < QMCKL_INVALID_EXIT_CODE);
-
-  qmckl_context_struct* ctx = (qmckl_context_struct*) context;
-  if (ctx == NULL) return QMCKL_FAILURE;
-
-  if (ctx->error != NULL) {
-    free(ctx->error);
-    ctx->error = NULL;
-  }
-
-  qmckl_error_struct* error = (qmckl_error_struct*) qmckl_malloc (context, sizeof(qmckl_error_struct));
-  error->exit_code = exit_code;
-  strcpy(error->function, function);
-  strcpy(error->message, message);
-
-  ctx->error = error;
-
-  return QMCKL_SUCCESS;
-}
-
      -
      -
      -
      2. - -
    2. TODO Test
      3. -
    -
    4. - -
  4. qmckl_context_set_error
    -
    -
    -
    qmckl_context
-qmckl_context_set_error(qmckl_context context, const qmckl_exit_code exit_code, const char* function, const char* message);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_context
-qmckl_context_set_error(qmckl_context context, const qmckl_exit_code exit_code, const char* function, const char* message)
-{
-  assert (context != 0);
-  assert (function != NULL);
-  assert (message != NULL);
-  assert (exit_code > 0);
-  assert (exit_code < QMCKL_INVALID_EXIT_CODE);
-
-  qmckl_context new_context = qmckl_context_copy(context);
-  if (new_context == 0) return context;
-
-  if (qmckl_context_update_error(new_context, exit_code,
-                                 function, message) != QMCKL_SUCCESS) {
-    return context;
-  }
-
-  return new_context;
-}
-
      -
      -
      -
      2. - -
    2. TODO Test
      3. -
    -
    5. - - -
  5. qmckl_context_check
    -
    -

    -Checks if the domain pointed by the pointer is a valid context. -Returns the input qmckl_context if the context is valid, 0 -otherwise. -

    - -
    -
    qmckl_context qmckl_context_check(const qmckl_context context) ;
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_context qmckl_context_check(const qmckl_context context) {
-
-  if (context == (qmckl_context) 0) return (qmckl_context) 0;
-
-  const qmckl_context_struct * ctx = (qmckl_context_struct*) context;
-
-  if (ctx->tag != VALID_TAG) return (qmckl_context) 0;
-
-  return context;
-}
-
      -
      -
      -
      2. -
    -
    6. - -
  6. qmckl_context_create
    -
    -

    -To create a new context, use qmckl_context_create(). -

    -
      -
    • On success, returns a pointer to a context using the qmckl_context type
      • -

    • -Returns 0 upon failure to allocate the internal data structure -

      - -
      -
      qmckl_context qmckl_context_create();
-
      -
      • -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_context qmckl_context_create() {
-
-  qmckl_context_struct* context =
-    (qmckl_context_struct*) qmckl_malloc ((qmckl_context) 0, sizeof(qmckl_context_struct));
-  if (context == NULL) {
-    return (qmckl_context) 0;
-  }
-
-  context->prev      = NULL;
-  context->ao_basis  = NULL;
-  context->precision = QMCKL_DEFAULT_PRECISION;
-  context->range     = QMCKL_DEFAULT_RANGE;
-  context->tag       = VALID_TAG;
-  context->error     = NULL;
-
-  return (qmckl_context) context;
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int64_t) function qmckl_context_create() bind(C)
-     use, intrinsic :: iso_c_binding
-   end function qmckl_context_create
-end interface
-
      -
      -
      -
      3. -
    -
    7. - -
  7. qmckl_context_copy
    -
    -

    -This function makes a shallow copy of the current context. -

    -
      -
    • Copying the 0-valued context returns 0
      • -
    • On success, returns a pointer to the new context using the qmckl_context type
      • -

    • -Returns 0 upon failure to allocate the internal data structure -for the new context -

      - -
      -
      qmckl_context qmckl_context_copy(const qmckl_context context);
-
      -
      • -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_context qmckl_context_copy(const qmckl_context context) {
-
-  const qmckl_context checked_context = qmckl_context_check(context);
-
-  if (checked_context == (qmckl_context) 0) {
-    return (qmckl_context) 0;
-  }
-
-  qmckl_context_struct* old_context = (qmckl_context_struct*) checked_context;
-
-  qmckl_context_struct* new_context =
-    (qmckl_context_struct*) qmckl_malloc (context, sizeof(qmckl_context_struct));
-  if (new_context == NULL) {
-    return (qmckl_context) 0;
-  }
-
-  new_context->prev      = old_context;
-  new_context->ao_basis  = old_context->ao_basis;
-  new_context->precision = old_context->precision;
-  new_context->range     = old_context->range;
-  new_context->tag       = VALID_TAG;
-  new_context->error     = old_context->error;
-
-  return (qmckl_context) new_context;
-}
-
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int64_t) function qmckl_context_copy(context) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-   end function qmckl_context_copy
-end interface
-
      -
      -
      -
      3. -
    -
    8. - -
  8. qmckl_context_previous
    -
    -

    -Returns the previous context -

    -
      -
    • On success, returns the ancestor of the current context
      • -
    • Returns 0 for the initial context
      • -

    • -Returns 0 for the 0-valued context -

      - -
      -
      qmckl_context qmckl_context_previous(const qmckl_context context);
-
      -
      • -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_context qmckl_context_previous(const qmckl_context context) {
-
-  const qmckl_context checked_context = qmckl_context_check(context);
-  if (checked_context == (qmckl_context) 0) {
-    return (qmckl_context) 0;
-  }
-
-  const qmckl_context_struct* ctx = (qmckl_context_struct*) checked_context;
-  return qmckl_context_check((qmckl_context) ctx->prev);
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int64_t) function qmckl_context_previous(context) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-   end function qmckl_context_previous
-end interface
-
      -
      -
      -
      3. -
    -
    9. - -
  9. qmckl_context_destroy
    -
    -

    -Destroys the current context, leaving the ancestors untouched. -

    -
      -
    • Succeeds if the current context is properly destroyed
      • -
    • Fails otherwise
      • -
    • Fails if the 0-valued context is given in argument
      • -
    • Fails if the the pointer is not a valid context
      • -
    - -
    -
    qmckl_exit_code qmckl_context_destroy(qmckl_context context);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_exit_code qmckl_context_destroy(const qmckl_context context) {
-
-  const qmckl_context checked_context = qmckl_context_check(context);
-  if (checked_context == (qmckl_context) 0) return QMCKL_FAILURE;
-
-  qmckl_context_struct* ctx = (qmckl_context_struct*) context;
-  if (ctx == NULL) return QMCKL_FAILURE;
-
-  ctx->tag = INVALID_TAG;
-  return qmckl_free(context,ctx);
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int32_t) function qmckl_context_destroy(context) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-   end function qmckl_context_destroy
-end interface
-
      -
      -
      -
      3. -
    -
    10. - -
  10. Basis set
    -
    -

    -For H2 with the following basis set, -

    - -
    -HYDROGEN
-S   5
-1         3.387000E+01           6.068000E-03
-2         5.095000E+00           4.530800E-02
-3         1.159000E+00           2.028220E-01
-4         3.258000E-01           5.039030E-01
-5         1.027000E-01           3.834210E-01
-S   1
-1         3.258000E-01           1.000000E+00
-S   1
-1         1.027000E-01           1.000000E+00
-P   1
-1         1.407000E+00           1.000000E+00
-P   1
-1         3.880000E-01           1.000000E+00
-D   1
-1         1.057000E+00           1.0000000
-
    - -

    -we have: -

    - -
    -type = 'G'
-shell_num = 12
-prim_num = 20
-SHELL_CENTER = [1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2]
-SHELL_ANG_MOM = ['S', 'S', 'S', 'P', 'P', 'D', 'S', 'S', 'S', 'P', 'P', 'D']
-SHELL_PRIM_NUM = [5, 1, 1, 1, 1, 1, 5, 1, 1, 1, 1, 1]
-prim_index = [1, 6, 7, 8, 9, 10, 11, 16, 17, 18, 19, 20]
-EXPONENT = [ 33.87, 5.095, 1.159, 0.3258, 0.1027, 0.3258, 0.1027,
-             1.407, 0.388, 1.057, 33.87, 5.095, 1.159, 0.3258, 0.1027,
-             0.3258, 0.1027, 1.407, 0.388, 1.057]
-COEFFICIENT = [ 0.006068, 0.045308, 0.202822, 0.503903, 0.383421,
-                1.0, 1.0, 1.0, 1.0, 1.0, 0.006068, 0.045308, 0.202822,
-                0.503903, 0.383421, 1.0, 1.0, 1.0, 1.0, 1.0]
-
    -
    -
    11. - -
  11. qmckl_context_update_ao_basis
    -
    -

    -Updates the data describing the AO basis set into the context. -

    - - - - --- -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    typeGaussian or Slater
    shell_numNumber of shells
    prim_numTotal number of primitives
    SHELL_CENTER(shell_num)Id of the nucleus on which the shell is centered
    SHELL_ANG_MOM(shell_num)Id of the nucleus on which the shell is centered
    SHELL_FACTOR(shell_num)Normalization factor for the shell
    SHELL_PRIM_NUM(shell_num)Number of primitives in the shell
    SHELL_PRIM_INDEX(shell_num)Address of the first primitive of the shelll in the EXPONENT array
    EXPONENT(prim_num)Array of exponents
    COEFFICIENT(prim_num)Array of coefficients
    - -
    -
    qmckl_exit_code
-qmckl_context_update_ao_basis(qmckl_context   context     , const char      type,
-                              const int64_t   shell_num   , const int64_t   prim_num,
-                              const int64_t * SHELL_CENTER, const int32_t * SHELL_ANG_MOM,
-                              const double  * SHELL_FACTOR, const int64_t * SHELL_PRIM_NUM,
-                              const int64_t * SHELL_PRIM_INDEX,
-                              const double  * EXPONENT    , const double  * COEFFICIENT);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_exit_code
-qmckl_context_update_ao_basis(qmckl_context   context     , const char      type,
-                              const int64_t   shell_num   , const int64_t   prim_num,
-                              const int64_t * SHELL_CENTER, const int32_t * SHELL_ANG_MOM,
-                              const double  * SHELL_FACTOR, const int64_t * SHELL_PRIM_NUM,
-                              const int64_t * SHELL_PRIM_INDEX,
-                              const double  * EXPONENT    , const double  * COEFFICIENT)
-{
-
-  int64_t i;
-
-  /* Check input */
-
-  if (type != 'G' && type != 'S') return QMCKL_FAILURE;
-  if (shell_num <= 0) return QMCKL_FAILURE;
-  if (prim_num  <= 0) return QMCKL_FAILURE;
-  if (prim_num  < shell_num) return QMCKL_FAILURE;
-
-  for (i=0 ; i<shell_num ; i++) {
-    if (SHELL_CENTER[i] <= 0) return QMCKL_FAILURE;
-    if (SHELL_PRIM_NUM[i] <= 0) return QMCKL_FAILURE;
-    if (SHELL_ANG_MOM[i] < 0) return QMCKL_FAILURE;
-    if (SHELL_PRIM_INDEX[i] < 0) return QMCKL_FAILURE;
-  }
-
-  for (i=0 ; i<prim_num ; i++) {
-    if (EXPONENT[i] <= 0) return QMCKL_FAILURE;
-  }
-
-  qmckl_context_struct* ctx = (qmckl_context_struct*) context;
-  if (ctx == NULL) return QMCKL_FAILURE;
-
-  qmckl_ao_basis_struct* basis = (qmckl_ao_basis_struct*) malloc (sizeof(qmckl_ao_basis_struct));
-  if (basis == NULL) return QMCKL_FAILURE;
-
-
-  /* Memory allocations */
-
-  basis->shell_center  = (int64_t*) malloc (shell_num * sizeof(int64_t));
-  if (basis->shell_center == NULL) {
-    qmckl_free(context, basis);
-    return QMCKL_FAILURE;
-  }
-
-  basis->shell_ang_mom = (int32_t*) malloc (shell_num * sizeof(int32_t));
-  if (basis->shell_ang_mom == NULL) {
-    qmckl_free(context, basis->shell_center);
-    qmckl_free(context, basis);
-    return QMCKL_FAILURE;
-  }
-
-  basis->shell_prim_num= (int64_t*) malloc (shell_num * sizeof(int64_t));
-  if (basis->shell_prim_num == NULL) {
-    qmckl_free(context, basis->shell_ang_mom);
-    qmckl_free(context, basis->shell_center);
-    qmckl_free(context, basis);
-    return QMCKL_FAILURE;
-  }
-
-  basis->shell_factor  = (double *) malloc (shell_num * sizeof(double ));
-  if (basis->shell_factor == NULL) {
-    qmckl_free(context, basis->shell_prim_num);
-    qmckl_free(context, basis->shell_ang_mom);
-    qmckl_free(context, basis->shell_center);
-    qmckl_free(context, basis);
-    return QMCKL_FAILURE;
-  }
-
-  basis->exponent      = (double *) malloc (prim_num  * sizeof(double ));
-  if (basis->exponent == NULL) {
-    qmckl_free(context, basis->shell_factor);
-    qmckl_free(context, basis->shell_prim_num);
-    qmckl_free(context, basis->shell_ang_mom);
-    qmckl_free(context, basis->shell_center);
-    qmckl_free(context, basis);
-    return QMCKL_FAILURE;
-  }
-
-  basis->coefficient   = (double *) malloc (prim_num  * sizeof(double ));
-  if (basis->coefficient == NULL) {
-    qmckl_free(context, basis->exponent);
-    qmckl_free(context, basis->shell_factor);
-    qmckl_free(context, basis->shell_prim_num);
-    qmckl_free(context, basis->shell_ang_mom);
-    qmckl_free(context, basis->shell_center);
-    qmckl_free(context, basis);
-    return QMCKL_FAILURE;
-  }
-
-
-  /* Assign data */
-
-  basis->type      = type;
-  basis->shell_num = shell_num;
-  basis->prim_num  = prim_num;
-
-  for (i=0 ; i<shell_num ; i++) {
-    basis->shell_center  [i] = SHELL_CENTER  [i];
-    basis->shell_ang_mom [i] = SHELL_ANG_MOM [i];
-    basis->shell_prim_num[i] = SHELL_PRIM_NUM[i];
-    basis->shell_factor  [i] = SHELL_FACTOR  [i];
-  }
-
-  for (i=0 ; i<prim_num ; i++) {
-    basis->exponent   [i] = EXPONENT[i];
-    basis->coefficient[i] = COEFFICIENT[i];
-  }
-
-  ctx->ao_basis = basis;
-  return QMCKL_SUCCESS;
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int32_t) function qmckl_context_update_ao_basis(context, &
-        typ, shell_num, prim_num, SHELL_CENTER, SHELL_ANG_MOM, SHELL_FACTOR, &
-        SHELL_PRIM_NUM, SHELL_PRIM_INDEX, EXPONENT, COEFFICIENT) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-     character(c_char)  , intent(in), value :: typ
-     integer (c_int64_t), intent(in), value :: shell_num
-     integer (c_int64_t), intent(in), value :: prim_num
-     integer (c_int64_t), intent(in)        :: SHELL_CENTER(shell_num)
-     integer (c_int32_t), intent(in)        :: SHELL_ANG_MOM(shell_num)
-     double precision   , intent(in)        :: SHELL_FACTOR(shell_num)
-     integer (c_int64_t), intent(in)        :: SHELL_PRIM_NUM(shell_num)
-     integer (c_int64_t), intent(in)        :: SHELL_PRIM_INDEX(shell_num)
-     double precision   , intent(in)        :: EXPONENT(prim_num)
-     double precision   , intent(in)        :: COEFFICIENT(prim_num)
-   end function qmckl_context_update_ao_basis
-end interface
-
      -
      -
      -
      3. - -
    3. TODO Test
      4. -
    -
    12. - -
  12. qmckl_context_set_ao_basis
    -
    -

    -Sets the data describing the AO basis set into the context. -

    - - - - --- -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    typeGaussian or Slater
    shell_numNumber of shells
    prim_numTotal number of primitives
    SHELL_CENTER(shell_num)Id of the nucleus on which the shell is centered
    SHELL_ANG_MOM(shell_num)Id of the nucleus on which the shell is centered
    SHELL_FACTOR(shell_num)Normalization factor for the shell
    SHELL_PRIM_NUM(shell_num)Number of primitives in the shell
    SHELL_PRIM_INDEX(shell_num)Address of the first primitive of the shelll in the EXPONENT array
    EXPONENT(prim_num)Array of exponents
    COEFFICIENT(prim_num)Array of coefficients
    - -
    -
    qmckl_context
-qmckl_context_set_ao_basis(const qmckl_context context     , const char      type,
-                           const int64_t       shell_num   , const int64_t   prim_num,
-                           const int64_t *     SHELL_CENTER, const int32_t * SHELL_ANG_MOM,
-                           const double  *     SHELL_FACTOR, const int64_t * SHELL_PRIM_NUM,
-                           const int64_t *     SHELL_PRIM_INDEX,
-                           const double  *     EXPONENT    , const double  * COEFFICIENT);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_context
-qmckl_context_set_ao_basis(const qmckl_context context     , const char      type,
-                           const int64_t       shell_num   , const int64_t   prim_num,
-                           const int64_t *     SHELL_CENTER, const int32_t * SHELL_ANG_MOM,
-                           const double  *     SHELL_FACTOR, const int64_t * SHELL_PRIM_NUM,
-                           const int64_t *     SHELL_PRIM_INDEX,
-                           const double  *     EXPONENT    , const double  * COEFFICIENT)
-{
-
-  qmckl_context new_context = qmckl_context_copy(context);
-  if (new_context == 0) return 0;
-
-  if (qmckl_context_update_ao_basis(new_context, type, shell_num, prim_num,
-                                    SHELL_CENTER, SHELL_ANG_MOM, SHELL_FACTOR,
-                                    SHELL_PRIM_NUM, SHELL_PRIM_INDEX, EXPONENT,
-                                    COEFFICIENT
-                                    ) == QMCKL_FAILURE)
-    return 0;
-
-  return new_context;
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int64_t) function qmckl_context_set_ao_basis(context, &
-        typ, shell_num, prim_num, SHELL_CENTER, SHELL_ANG_MOM, SHELL_FACTOR, &
-        SHELL_PRIM_NUM, SHELL_PRIM_INDEX, EXPONENT, COEFFICIENT) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-     character(c_char)  , intent(in), value :: typ
-     integer (c_int64_t), intent(in), value :: shell_num
-     integer (c_int64_t), intent(in), value :: prim_num
-     integer (c_int64_t), intent(in)        :: SHELL_CENTER(shell_num)
-     integer (c_int32_t), intent(in)        :: SHELL_ANG_MOM(shell_num)
-     double precision   , intent(in)        :: SHELL_FACTOR(shell_num)
-     integer (c_int64_t), intent(in)        :: SHELL_PRIM_NUM(shell_num)
-     integer (c_int64_t), intent(in)        :: SHELL_PRIM_INDEX(shell_num)
-     double precision   , intent(in)        :: EXPONENT(prim_num)
-     double precision   , intent(in)        :: COEFFICIENT(prim_num)
-   end function qmckl_context_set_ao_basis
-end interface
-
      -
      -
      -
      3. - -
    3. TODO Test
      4. -
    -
    13. - -
  13. Precision
    -
    -

    -The following functions set and get the expected required -precision and range. precision should be an integer between 2 -and 53, and range should be an integer between 2 and 11. -

    - -

    -The setter functions functions return a new context as a 64-bit -integer. The getter functions return the value, as a 32-bit -integer. The update functions return QMCKL_SUCCESS or -QMCKL_FAILURE. -

    -
    -
    14. - -
  14. qmckl_context_update_precision
    -
    -

    -Modifies the parameter for the numerical precision in a given context. -

    -
    -
    qmckl_exit_code qmckl_context_update_precision(const qmckl_context context, const int precision);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_exit_code qmckl_context_update_precision(const qmckl_context context, const int precision) {
-
-  if (precision <  2) return QMCKL_FAILURE;
-  if (precision > 53) return QMCKL_FAILURE;
-
-  qmckl_context_struct* ctx = (qmckl_context_struct*) context;
-  if (ctx == NULL) return QMCKL_FAILURE;
-
-  ctx->precision = precision;
-  return QMCKL_SUCCESS;
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int32_t) function qmckl_context_update_precision(context, precision) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-     integer (c_int32_t), intent(in), value :: precision
-   end function qmckl_context_update_precision
-end interface
-
      -
      -
      -
      3. -
    -
    15. -
  15. qmckl_context_update_range
    -
    -

    -Modifies the parameter for the numerical range in a given context. -

    -
    -
    qmckl_exit_code qmckl_context_update_range(const qmckl_context context, const int range);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_exit_code qmckl_context_update_range(const qmckl_context context, const int range) {
-
-  if (range <  2) return QMCKL_FAILURE;
-  if (range > 11) return QMCKL_FAILURE;
-
-  qmckl_context_struct* ctx = (qmckl_context_struct*) context;
-  if (ctx == NULL) return QMCKL_FAILURE;
-
-  ctx->range = range;
-  return QMCKL_SUCCESS;
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int32_t) function qmckl_context_update_range(context, range) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-     integer (c_int32_t), intent(in), value :: range
-   end function qmckl_context_update_range
-end interface
-
      -
      -
      -
      3. -
    -
    16. -
  16. qmckl_context_set_precision
    -
    -

    -Returns a copy of the context with a different precision parameter. -

    -
    -
    qmckl_context qmckl_context_set_precision(const qmckl_context context, const int precision);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_context qmckl_context_set_precision(const qmckl_context context, const int precision) {
-  qmckl_context new_context = qmckl_context_copy(context);
-  if (new_context == 0) return 0;
-
-  if (qmckl_context_update_precision(new_context, precision) == QMCKL_FAILURE) return 0;
-
-  return new_context;
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int64_t) function qmckl_context_set_precision(context, precision) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-     integer (c_int32_t), intent(in), value :: precision
-   end function qmckl_context_set_precision
-end interface
-
      -
      -
      -
      3. -
    -
    17. -
  17. qmckl_context_set_range
    -
    -

    -Returns a copy of the context with a different precision parameter. -

    -
    -
    qmckl_context qmckl_context_set_range(const qmckl_context context, const int range);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      qmckl_context qmckl_context_set_range(const qmckl_context context, const int range) {
-  qmckl_context new_context = qmckl_context_copy(context);
-  if (new_context == 0) return 0;
-
-  if (qmckl_context_update_range(new_context, range) == QMCKL_FAILURE) return 0;
-
-  return new_context;
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int64_t) function qmckl_context_set_range(context, range) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-     integer (c_int32_t), intent(in), value :: range
-   end function qmckl_context_set_range
-end interface
-
      -
      -
      -
      3. -
    -
    18. - -
  18. qmckl_context_get_precision
    -
    -

    -Returns the value of the numerical precision in the context -

    -
    -
    int32_t qmckl_context_get_precision(const qmckl_context context);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      int qmckl_context_get_precision(const qmckl_context context) {
-  const qmckl_context_struct* ctx = (qmckl_context_struct*) context;
-  return ctx->precision;
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int32_t) function qmckl_context_get_precision(context) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-   end function qmckl_context_get_precision
-end interface
-
      -
      -
      -
      3. -
    -
    19. -
  19. qmckl_context_get_range
    -
    -

    -Returns the value of the numerical range in the context -

    -
    -
    int32_t qmckl_context_get_range(const qmckl_context context);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      int qmckl_context_get_range(const qmckl_context context) {
-  const qmckl_context_struct* ctx = (qmckl_context_struct*) context;
-  return ctx->range;
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   integer (c_int32_t) function qmckl_context_get_range(context) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-   end function qmckl_context_get_range
-end interface
-
      -
      -
      -
      3. -
    -
    20. - -
  20. qmckl_context_get_epsilon
    -
    -

    -Returns \(\epsilon = 2^{1-n}\) where n is the precision -

    -
    -
    double qmckl_context_get_epsilon(const qmckl_context context);
-
    -
    -
    - -
      -
    1. Source
      -
      -
      -
      double qmckl_context_get_epsilon(const qmckl_context context) {
-  const qmckl_context_struct* ctx = (qmckl_context_struct*) context;
-  return pow(2.0,(double) 1-ctx->precision);
-}
-
      -
      -
      -
      2. - -
    2. Fortran interface
      -
      -
      -
      interface
-   real (c_double) function qmckl_context_get_epsilon(context) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-   end function qmckl_context_get_epsilon
-end interface
-
      -
      -
      -
      3. -
    -
    21. -
-
-
-
-

2.3 Error handling

-
-

-This file is written in C because it is more natural to express the -error handling in C than in Fortran. -

- -

-2 files are produced: -

-
    -
  • a source file : qmckl_error.c
    • -
  • a test file : test_qmckl_error.c
    • -
-
- -
-

2.3.1 Error handling

-
-

-The library should never make the calling programs abort, nor -perform any input/output operations. This decision has to be taken -by the developer of the code calling the library. -

- -

-All the functions return with an exit code, defined as -

-
-
typedef int32_t qmckl_exit_code;
-
-
- -

-The exit code returns the completion status of the function to the -calling program. When a function call completed successfully, the -QMCKL_SUCCESS exit code is returned. If one of the functions of -the library fails to complete the requested task, an appropriate -error code is returned to the program. -

- -

-Here is the complete list of exit codes. -

- - - - --- -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
QMCKL_SUCCESS0
QMCKL_INVALID_ARG_11
QMCKL_INVALID_ARG_22
QMCKL_INVALID_ARG_33
QMCKL_INVALID_ARG_44
QMCKL_INVALID_ARG_55
QMCKL_INVALID_ARG_66
QMCKL_INVALID_ARG_77
QMCKL_INVALID_ARG_88
QMCKL_INVALID_ARG_99
QMCKL_INVALID_ARG_1010
QMCKL_NULL_CONTEXT101
QMCKL_FAILURE102
QMCKL_ERRNO103
QMCKL_INVALID_EXIT_CODE104
- -
-
""" This script generates the C and Fortran constants for the error
-    codes from the org-mode table.
-"""
-
-result = [ "#+BEGIN_SRC C :comments org :tangle qmckl.h" ] 
-for (text, code) in table:
-    text=text.replace("~","")
-    result += [ f"#define  {text:30s} {code:d}" ]
-result += [ "#+END_SRC" ]
-    
-result += [ "" ] 
-
-result += [ "#+BEGIN_SRC f90 :comments org :tangle qmckl_f.f90" ]
-for (text, code) in table:
-    text=text.replace("~","")
-    result += [ f"   integer, parameter :: {text:30s} = {code:d}" ]
-result += [ "#+END_SRC" ]
-
-return '\n'.join(result)
-
-
-
- -
-
#define  QMCKL_SUCCESS                  0
-#define  QMCKL_INVALID_ARG_1            1
-#define  QMCKL_INVALID_ARG_2            2
-#define  QMCKL_INVALID_ARG_3            3
-#define  QMCKL_INVALID_ARG_4            4
-#define  QMCKL_INVALID_ARG_5            5
-#define  QMCKL_INVALID_ARG_6            6
-#define  QMCKL_INVALID_ARG_7            7
-#define  QMCKL_INVALID_ARG_8            8
-#define  QMCKL_INVALID_ARG_9            9
-#define  QMCKL_INVALID_ARG_10           10
-#define  QMCKL_NULL_CONTEXT             101
-#define  QMCKL_FAILURE                  102
-#define  QMCKL_ERRNO                    103
-#define  QMCKL_INVALID_EXIT_CODE        104
-
-
- -
-
integer, parameter :: QMCKL_SUCCESS                  = 0
-integer, parameter :: QMCKL_INVALID_ARG_1            = 1
-integer, parameter :: QMCKL_INVALID_ARG_2            = 2
-integer, parameter :: QMCKL_INVALID_ARG_3            = 3
-integer, parameter :: QMCKL_INVALID_ARG_4            = 4
-integer, parameter :: QMCKL_INVALID_ARG_5            = 5
-integer, parameter :: QMCKL_INVALID_ARG_6            = 6
-integer, parameter :: QMCKL_INVALID_ARG_7            = 7
-integer, parameter :: QMCKL_INVALID_ARG_8            = 8
-integer, parameter :: QMCKL_INVALID_ARG_9            = 9
-integer, parameter :: QMCKL_INVALID_ARG_10           = 10
-integer, parameter :: QMCKL_NULL_CONTEXT             = 101
-integer, parameter :: QMCKL_FAILURE                  = 102
-integer, parameter :: QMCKL_ERRNO                    = 103
-integer, parameter :: QMCKL_INVALID_EXIT_CODE        = 104
-
-
- -

-To make a function fail, the qmckl_failwith function should be -called, such that information about the failure is stored in -the context. The desired exit code is given as an argument, as -well as the name of the function and an error message. The return -code of the function is the desired return code. -

- -
-
qmckl_exit_code qmckl_failwith(qmckl_context context,
-                               const qmckl_exit_code exit_code,
-                               const char* function,
-                               const char* message) ;
-
-
- -
-
qmckl_exit_code qmckl_failwith(qmckl_context context,
-                               const qmckl_exit_code exit_code,
-                               const char* function,
-                               const char* message) {
-  if (context == 0) return QMCKL_NULL_CONTEXT;
-  assert (exit_code > 0);
-  assert (exit_code < QMCKL_INVALID_EXIT_CODE);
-  assert (function != NULL);
-  assert (message  != NULL);
-  assert (strlen(function) < QMCKL_MAX_FUN_LEN);
-  assert (strlen(message)  < QMCKL_MAX_MSG_LEN);
-
-  context = qmckl_context_set_error(context, exit_code, function, message);
-  return exit_code;
-}
-
-
-
- -

-For example, this function can be used as -

-
-
if (x < 0) {
-  return qmckl_failwith(context,
-                        QMCKL_INVALID_ARG_2,
-                        "qmckl_function", 
-                        "Expected x >= 0");
- }
-
-
-
-
- -
-

2.3.2 Multi-precision related constants

-
-

-Controlling numerical precision enables optimizations. Here, the -default parameters determining the target numerical precision and -range are defined. -

- -
-
#define QMCKL_DEFAULT_PRECISION 53
-#define QMCKL_DEFAULT_RANGE     11
-
-
- -
-
integer, parameter :: QMCKL_DEFAULT_PRECISION = 53
-integer, parameter :: QMCKL_DEFAULT_RANGE = 11
-
-
-
-
-
-
-

2.4 Memory management

-
-

-We override the allocation functions to enable the possibility of -optimized libraries to fine-tune the memory allocation. -

- -

-2 files are produced: -

-
    -
  • a source file : qmckl_memory.c
    • -
  • a test file : test_qmckl_memory.c
    • -
-
- -
-

2.4.1 qmckl_malloc

-
-

-Memory allocation function, letting the library choose how the -memory will be allocated, and a pointer is returned to the user. -The context is passed to let the library store data related to the -allocation inside the context. -

- -
-
void* qmckl_malloc(const qmckl_context ctx, const size_t size);
-
-
- -
-
interface
-   type (c_ptr) function qmckl_malloc (context, size) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-     integer (c_int64_t), intent(in), value :: size
-   end function qmckl_malloc
-end interface
-
-
-
- -
    -
  1. Source
    -
    -
    -
    void* qmckl_malloc(const qmckl_context ctx, const size_t size) {
-  if (ctx == (qmckl_context) 0) {}; /* Avoid unused argument warning */
-  void * result = malloc( (size_t) size );
-  assert (result != NULL) ;
-  return result;
-}
-
-
    -
    -
    -
    2. -
-
- -
-

2.4.2 qmckl_free

-
-

-The context is passed, in case some important information has been -stored related to memory allocation and needs to be updated. -

- -
-
qmckl_exit_code qmckl_free(qmckl_context context, void *ptr);
-
-
- -
-
interface
-   integer (c_int32_t) function qmckl_free (context, ptr) bind(C)
-     use, intrinsic :: iso_c_binding
-     integer (c_int64_t), intent(in), value :: context
-     type (c_ptr), intent(in), value :: ptr
-   end function qmckl_free
-end interface
-
-
-
- -
    -
  1. Source
    -
    -
    -
    qmckl_exit_code qmckl_free(qmckl_context context, void *ptr) {
-  if (context == 0) return QMCKL_INVALID_ARG_1;
-  if (ptr == NULL)  return QMCKL_INVALID_ARG_2;
-  free(ptr);
-  return QMCKL_SUCCESS;
-}
-
    -
    -
    -
    2. -
-
-
-
-

2.5 Computation of distances

-
-

-Function for the computation of distances between particles. -

- -

-3 files are produced: -

-
    -
  • a source file : qmckl_distance.f90
    • -
  • a C test file : test_qmckl_distance.c
    • -
  • a Fortran test file : test_qmckl_distance_f.f90
    • -
-
- -
-

2.5.1 Squared distance

-
-
-
    -
  1. qmckl_distance_sq
    -
    -

    -Computes the matrix of the squared distances between all pairs of -points in two sets, one point within each set: -\[ - C_{ij} = \sum_{k=1}^3 (A_{k,i}-B_{k,j})^2 - \] -

    -
    - -
      -
    1. Arguments
      -
      - - - --- -- -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      contextinputGlobal state
      transainputArray A is N: Normal, T: Transposed
      transbinputArray B is N: Normal, T: Transposed
      minputNumber of points in the first set
      ninputNumber of points in the second set
      A(lda,3)inputArray containing the \(m \times 3\) matrix \(A\)
      ldainputLeading dimension of array A
      B(ldb,3)inputArray containing the \(n \times 3\) matrix \(B\)
      ldbinputLeading dimension of array B
      C(ldc,n)outputArray containing the \(m \times n\) matrix \(C\)
      ldcinputLeading dimension of array C
      -
      -
      2. - -
    2. Requirements
      -
      -
        -
      • context is not 0
        • -
      • m > 0
        • -
      • n > 0
        • -
      • lda >= 3 if transa is N
        • -
      • lda >= m if transa is T
        • -
      • ldb >= 3 if transb is N
        • -
      • ldb >= n if transb is T
        • -
      • ldc >= m
        • -
      • A is allocated with at least \(3 \times m \times 8\) bytes
        • -
      • B is allocated with at least \(3 \times n \times 8\) bytes
        • -
      • C is allocated with at least \(m \times n \times 8\) bytes
        • -
      -
      -
      3. - -
    3. Performance
      -
      -

      -This function might be more efficient when A and B are -transposed. -

      - -
      -
      qmckl_exit_code qmckl_distance_sq(const qmckl_context context,
-                                  const char transa, const char transb,
-                                  const int64_t m, const int64_t n,
-                                  const double *A, const int64_t lda,
-                                  const double *B, const int64_t ldb,
-                                  const double *C, const int64_t ldc);
-
      -
      -
      -
      4. - -
    4. Source
      -
      -
      -
      integer function qmckl_distance_sq_f(context, transa, transb, m, n, A, LDA, B, LDB, C, LDC) result(info)
-  implicit none
-  integer*8  , intent(in)  :: context
-  character  , intent(in)  :: transa, transb
-  integer*8  , intent(in)  :: m, n
-  integer*8  , intent(in)  :: lda
-  real*8     , intent(in)  :: A(lda,*)
-  integer*8  , intent(in)  :: ldb
-  real*8     , intent(in)  :: B(ldb,*)
-  integer*8  , intent(in)  :: ldc
-  real*8     , intent(out) :: C(ldc,*)
-
-  integer*8 :: i,j
-  real*8    :: x, y, z
-  integer   :: transab
-
-  info = 0
-
-  if (context == 0_8) then
-     info = -1
-     return
-  endif
-
-  if (m <= 0_8) then
-     info = -2
-     return
-  endif
-
-  if (n <= 0_8) then
-     info = -3
-     return
-  endif
-
-  if (transa == 'N' .or. transa == 'n') then
-     transab = 0
-  else if (transa == 'T' .or. transa == 't') then
-     transab = 1
-  else
-     transab = -100
-  endif
-
-  if (transb == 'N' .or. transb == 'n') then
-     continue
-  else if (transa == 'T' .or. transa == 't') then
-     transab = transab + 2
-  else
-     transab = -100
-  endif
-
-  if (transab < 0) then
-     info = -4
-     return 
-  endif
-
-  if (iand(transab,1) == 0 .and. LDA < 3) then
-     info = -5
-     return
-  endif
-
-  if (iand(transab,1) == 1 .and. LDA < m) then
-     info = -6
-     return
-  endif
-
-  if (iand(transab,2) == 0 .and. LDA < 3) then
-     info = -6
-     return
-  endif
-
-  if (iand(transab,2) == 2 .and. LDA < m) then
-     info = -7
-     return
-  endif
-
-
-  select case (transab)
-
-  case(0)
-
-     do j=1,n
-        do i=1,m
-           x = A(1,i) - B(1,j)
-           y = A(2,i) - B(2,j)
-           z = A(3,i) - B(3,j)
-           C(i,j) = x*x + y*y + z*z
-        end do
-     end do
-
-  case(1)
-
-     do j=1,n
-        do i=1,m
-           x = A(i,1) - B(1,j)
-           y = A(i,2) - B(2,j)
-           z = A(i,3) - B(3,j)
-           C(i,j) = x*x + y*y + z*z
-        end do
-     end do
-
-  case(2)
-
-     do j=1,n
-        do i=1,m
-           x = A(1,i) - B(j,1)
-           y = A(2,i) - B(j,2)
-           z = A(3,i) - B(j,3)
-           C(i,j) = x*x + y*y + z*z
-        end do
-     end do
-
-  case(3)
-
-     do j=1,n
-        do i=1,m
-           x = A(i,1) - B(j,1)
-           y = A(i,2) - B(j,2)
-           z = A(i,3) - B(j,3)
-           C(i,j) = x*x + y*y + z*z
-        end do
-     end do
-
-  end select
-
-end function qmckl_distance_sq_f
-
      -
      -
      -
      5. -
    -
    2. -
-
-
-
-

2.6 Atomic Orbitals

-
-

-This files contains all the routines for the computation of the -values, gradients and Laplacian of the atomic basis functions. -

- -

-3 files are produced: -

-
    -
  • a source file : qmckl_ao.f90
    • -
  • a C test file : test_qmckl_ao.c
    • -
  • a Fortran test file : test_qmckl_ao_f.f90
    • -
-
- -
-

2.6.1 Polynomials

-
-

-\[ - P_l(\mathbf{r},\mathbf{R}_i) = (x-X_i)^a (y-Y_i)^b (z-Z_i)^c - \] -

- -\begin{eqnarray*} -\frac{\partial }{\partial x} P_l\left(\mathbf{r},\mathbf{R}_i \right) & - = & a (x-X_i)^{a-1} (y-Y_i)^b (z-Z_i)^c \\ -\frac{\partial }{\partial y} P_l\left(\mathbf{r},\mathbf{R}_i \right) & - = & b (x-X_i)^a (y-Y_i)^{b-1} (z-Z_i)^c \\ -\frac{\partial }{\partial z} P_l\left(\mathbf{r},\mathbf{R}_i \right) & - = & c (x-X_i)^a (y-Y_i)^b (z-Z_i)^{c-1} \\ -\end{eqnarray*} - - -\begin{eqnarray*} -\left( \frac{\partial }{\partial x^2} + - \frac{\partial }{\partial y^2} + - \frac{\partial }{\partial z^2} \right) P_l - \left(\mathbf{r},\mathbf{R}_i \right) & = & - a(a-1) (x-X_i)^{a-2} (y-Y_i)^b (z-Z_i)^c + \\ - && b(b-1) (x-X_i)^a (y-Y_i)^{b-1} (z-Z_i)^c + \\ - && c(c-1) (x-X_i)^a (y-Y_i)^b (z-Z_i)^{c-1} -\end{eqnarray*} -
- -
    -
  1. qmckl_ao_power
    -
    -

    -Computes all the powers of the n input data up to the given -maximum value given in input for each of the \(n\) points: -

    - -

    -\[ P_{ij} = X_j^i \] -

    -
    - -
      -
    1. Arguments
      -
      - - - --- -- -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      contextinputGlobal state
      ninputNumber of values
      X(n)inputArray containing the input values
      LMAX(n)inputArray containing the maximum power for each value
      P(LDP,n)outputArray containing all the powers of X
      LDPinputLeading dimension of array P
      -
      -
      2. - -
    2. Requirements
      -
      -
        -
      • context is not 0
        • -
      • n > 0
        • -
      • X is allocated with at least \(n \times 8\) bytes
        • -
      • LMAX is allocated with at least \(n \times 4\) bytes
        • -
      • P is allocated with at least \(n \times \max_i \text{LMAX}_i \times 8\) bytes
        • -
      • LDP >= \(\max_i\) LMAX[i]
        • -
      -
      -
      3. - -
    3. Header
      -
      -
      -
      qmckl_exit_code qmckl_ao_power(const qmckl_context context,
-                const int64_t n, 
-                const double *X, const int32_t *LMAX,
-                const double *P, const int64_t LDP);
-
      -
      -
      -
      4. - -
    4. Source
      -
      -
      -
      integer function qmckl_ao_power_f(context, n, X, LMAX, P, ldp) result(info)
-  implicit none
-  integer*8 , intent(in)  :: context
-  integer*8 , intent(in)  :: n
-  real*8    , intent(in)  :: X(n)
-  integer   , intent(in)  :: LMAX(n)
-  real*8    , intent(out) :: P(ldp,n)
-  integer*8 , intent(in)  :: ldp
-
-  integer*8  :: i,j
-
-  info = 0
-
-  if (context == 0_8) then
-     info = -1
-     return
-  endif
-
-  if (LDP < MAXVAL(LMAX)) then
-     info = -2
-     return
-  endif
-
-  do j=1,n
-    P(1,j) = X(j)
-    do i=2,LMAX(j)
-       P(i,j) = P(i-1,j) * X(j) 
-    end do
-  end do
-
-end function qmckl_ao_power_f
-
      -
      -
      -
      5. -
    -
    2. - - -
  2. qmckl_ao_polynomial_vgl
    -
    -

    -Computes the values, gradients and Laplacians at a given point of -all polynomials with an angular momentum up to lmax. -

    -
    - -
      -
    1. Arguments
      -
      - - - --- -- -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      contextinputGlobal state
      X(3)inputArray containing the coordinates of the points
      R(3)inputArray containing the x,y,z coordinates of the center
      lmaxinputMaximum angular momentum
      noutputNumber of computed polynomials
      L(ldl,n)outputContains a,b,c for all n results
      ldlinputLeading dimension of L
      VGL(ldv,n)outputValue, gradients and Laplacian of the polynomials
      ldvinputLeading dimension of array VGL
      -
      -
      2. - -
    2. Requirements
      -
      -
        -
      • context is not 0
        • -
      • n > 0
        • -
      • lmax >= 0
        • -
      • ldl >= 3
        • -
      • ldv >= 5
        • -
      • X is allocated with at least \(3 \times 8\) bytes
        • -
      • R is allocated with at least \(3 \times 8\) bytes
        • -
      • n >= (lmax+1)(lmax+2)(lmax+3)/6
        • -
      • L is allocated with at least \(3 \times n \times 4\) bytes
        • -
      • VGL is allocated with at least \(5 \times n \times 8\) bytes
        • -
      • On output, n should be equal to (lmax+1)(lmax+2)(lmax+3)/6
        • -
      • On output, the powers are given in the following order (l=a+b+c): -
          -
        • Increase values of l
          • -
        • Within a given value of l, alphabetical order of the -string made by a*"x" + b*"y" + c*"z" (in Python notation). -For example, with a=0, b=2 and c=1 the string is "yyz"
          • -
        • -
      -
      -
      3. - -
    3. Error codes
      -
      - - - --- -- - - - - - - - - - - - - - - - - - - - - - -
      -1Null context
      -2Inconsistent ldl
      -3Inconsistent ldv
      -4Inconsistent lmax
      -
      -
      4. - -
    4. Header
      -
      -
      -
      qmckl_exit_code qmckl_ao_polynomial_vgl(const qmckl_context context,
-                const double *X, const double *R,
-                const int32_t lmax, const int64_t *n,
-                const int32_t *L,   const int64_t ldl,
-                const double *VGL,  const int64_t ldv);
-
      -
      -
      -
      5. - -
    5. Source
      -
      -
      -
      integer function qmckl_ao_polynomial_vgl_f(context, X, R, lmax, n, L, ldl, VGL, ldv) result(info)
-  implicit none
-  integer*8 , intent(in)  :: context
-  real*8    , intent(in)  :: X(3), R(3)
-  integer   , intent(in)  :: lmax
-  integer*8 , intent(out) :: n
-  integer   , intent(out) :: L(ldl,(lmax+1)*(lmax+2)*(lmax+3)/6)
-  integer*8 , intent(in)  :: ldl
-  real*8    , intent(out) :: VGL(ldv,(lmax+1)*(lmax+2)*(lmax+3)/6)
-  integer*8 , intent(in)  :: ldv
-
-  integer*8         :: i,j
-  integer           :: a,b,c,d
-  real*8            :: Y(3)
-  integer           :: lmax_array(3)
-  real*8            :: pows(-2:lmax,3)
-  integer, external :: qmckl_ao_power_f
-  double precision  :: xy, yz, xz
-  double precision  :: da, db, dc, dd
-
-  info = 0
-
-  if (context == 0_8) then
-     info = -1
-     return
-  endif
-
-  if (ldl < 3) then
-     info = -2
-     return
-  endif
-
-  if (ldv < 5) then
-     info = -3
-     return
-  endif
-
-  if (lmax <= 0) then
-     info = -4
-     return
-  endif
-
-
-  do i=1,3
-     Y(i) = X(i) - R(i)
-  end do
-
-  lmax_array(1:3) = lmax
-  if (lmax == 0) then
-    VGL(1,1) = 1.d0
-    vgL(2:5,1) = 0.d0
-    l(1:3,1) = 0
-    n=1
-  else if (lmax > 0) then
-    pows(-2:0,1:3) = 1.d0
-    do i=1,lmax
-        pows(i,1) = pows(i-1,1) * Y(1) 
-        pows(i,2) = pows(i-1,2) * Y(2) 
-        pows(i,3) = pows(i-1,3) * Y(3) 
-    end do
-
-    VGL(1:5,1:4) = 0.d0
-    l(1:3,1:4) = 0
-
-    VGL(1,1) = 1.d0
-    vgl(1:5,2:4) = 0.d0
-
-    l(1,2) = 1
-    vgl(1,2) = pows(1,1)
-    vgL(2,2) = 1.d0
-
-    l(2,3) = 1
-    vgl(1,3) = pows(1,2)
-    vgL(3,3) = 1.d0
-
-    l(3,4) = 1
-    vgl(1,4) = pows(1,3)
-    vgL(4,4) = 1.d0
-
-    n=4
-  endif
-
-  ! l>=2
-  dd = 2.d0
-  do d=2,lmax
-     da = dd
-     do a=d,0,-1
-        db = dd-da
-        do b=d-a,0,-1
-           c  = d  - a  - b
-           dc = dd - da - db
-           n = n+1
-
-           l(1,n) = a
-           l(2,n) = b
-           l(3,n) = c
-
-           xy = pows(a,1) * pows(b,2)
-           yz = pows(b,2) * pows(c,3)
-           xz = pows(a,1) * pows(c,3)
-
-           vgl(1,n) = xy * pows(c,3)
-
-           xy = dc * xy
-           xz = db * xz
-           yz = da * yz
-
-           vgl(2,n) = pows(a-1,1) * yz
-           vgl(3,n) = pows(b-1,2) * xz
-           vgl(4,n) = pows(c-1,3) * xy
-
-           vgl(5,n) = &
-                (da-1.d0) * pows(a-2,1) * yz + &
-                (db-1.d0) * pows(b-2,2) * xz + &
-                (dc-1.d0) * pows(c-2,3) * xy
-
-           db = db - 1.d0
-        end do
-        da = da - 1.d0
-     end do
-     dd = dd + 1.d0
-  end do
-
-  info = 0
-
-end function qmckl_ao_polynomial_vgl_f
-
      -
      -
      -
      6. -
    -
    3. -
-
- -
-

2.6.2 Gaussian basis functions

-
-
-
    -
  1. qmckl_ao_gaussian_vgl
    -
    -

    -Computes the values, gradients and Laplacians at a given point of -n Gaussian functions centered at the same point: +implementation of the library, we focus on the definition of the API +and the tests, and on a pedagogical presentation of the +algorithms. We expect the HPC experts to use this repository as a +reference for re-writing optimized libraries.

    -\[ v_i = \exp(-a_i |X-R|^2) \] -\[ \nabla_x v_i = -2 a_i (X_x - R_x) v_i \] -\[ \nabla_y v_i = -2 a_i (X_y - R_y) v_i \] -\[ \nabla_z v_i = -2 a_i (X_z - R_z) v_i \] -\[ \Delta v_i = a_i (4 |X-R|^2 a_i - 6) v_i \] +The source code of the library is available at +https://github.com/trex-coe/qmckl +and bug reports should be submitted at +https://github.com/trex-coe/qmckl/issues.

    -
    - -
      -
    1. Arguments
      -
      - - - --- -- -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      contextinputGlobal state
      X(3)inputArray containing the coordinates of the points
      R(3)inputArray containing the x,y,z coordinates of the center
      ninputNumber of computed gaussians
      A(n)inputExponents of the Gaussians
      VGL(ldv,5)outputValue, gradients and Laplacian of the Gaussians
      ldvinputLeading dimension of array VGL
      -
      -
      2. - -
    2. Requirements
      -
      -
        -
      • context is not 0
        • -
      • n > 0
        • -
      • ldv >= 5
        • -
      • A(i) > 0 for all i
        • -
      • X is allocated with at least \(3 \times 8\) bytes
        • -
      • R is allocated with at least \(3 \times 8\) bytes
        • -
      • A is allocated with at least \(n \times 8\) bytes
        • -
      • VGL is allocated with at least \(n \times 5 \times 8\) bytes
        • -
      -
      -
      3. - -
    3. Header
      -
      -
      -
      qmckl_exit_code qmckl_ao_gaussian_vgl(const qmckl_context context,
-                const double *X, const double *R,
-                const int64_t *n, const int64_t *A,
-                const double *VGL,  const int64_t ldv);
-
      -
      -
      -
      4. - -
    4. Source
      -
      -
      -
      integer function qmckl_ao_gaussian_vgl_f(context, X, R, n, A, VGL, ldv) result(info)
-  implicit none
-  integer*8 , intent(in)  :: context
-  real*8    , intent(in)  :: X(3), R(3)
-  integer*8 , intent(in)  :: n
-  real*8    , intent(in)  :: A(n)
-  real*8    , intent(out) :: VGL(ldv,5)
-  integer*8 , intent(in)  :: ldv
-
-  integer*8         :: i,j
-  real*8            :: Y(3), r2, t, u, v
-
-  info = 0
-
-  if (context == 0_8) then
-     info = -1
-     return
-  endif
-
-  if (n <= 0) then
-     info = -2
-     return
-  endif
-
-  if (ldv < n) then
-     info = -3
-     return
-  endif
-
-
-  do i=1,3
-     Y(i) = X(i) - R(i)
-  end do
-  r2 = Y(1)*Y(1) + Y(2)*Y(2) + Y(3)*Y(3)
-
-  do i=1,n
-     VGL(i,1) = dexp(-A(i) * r2)
-  end do
-
-  do i=1,n
-     VGL(i,5) = A(i) * VGL(i,1)
-  end do
-
-  t = -2.d0 * ( X(1) - R(1) )
-  u = -2.d0 * ( X(2) - R(2) )
-  v = -2.d0 * ( X(3) - R(3) )
-
-  do i=1,n
-     VGL(i,2) = t * VGL(i,5)
-     VGL(i,3) = u * VGL(i,5)
-     VGL(i,4) = v * VGL(i,5)
-  end do
-
-  t = 4.d0 * r2
-  do i=1,n
-     VGL(i,5) = (t * A(i) - 6.d0) *  VGL(i,5)
-  end do
-
-end function qmckl_ao_gaussian_vgl_f
-
      -
      -
      -
      5. -
    -
    2. -
-
- -
-

2.6.3 TODO Slater basis functions

-
-
-
+
-
-

3 Acknowledgments

-

-euflag.jpg -TREX: Targeting Real Chemical Accuracy at the Exascale project has received funding from the European Union’s Horizon 2020 - Research and Innovation program - under grant agreement no. 952165. The content of this document does not represent the opinion of the European Union, and the European Union is not responsible for any use that might be made of such content. +euflag.jpg TREX: Targeting Real Chemical Accuracy at the Exascale project has received funding from the European Union’s Horizon 2020 - Research and Innovation program - under grant agreement no. 952165. The content of this document does not represent the opinion of the European Union, and the European Union is not responsible for any use that might be made of such content.

-
-
-

Created: 2021-03-06 Sat 22:28

+

Author: TREX CoE

+

Created: 2021-03-19 Fri 12:49

Validate

diff --git a/qmckl.css b/qmckl.css new file mode 100644 index 0000000..cfc7a63 --- /dev/null +++ b/qmckl.css @@ -0,0 +1,972 @@ +/* Adapted from worg.css */ + +@import url(https://fonts.googleapis.com/css?family=Droid+Sans|Droid+Sans+Mono|Droid+Serif); + +@media all +{ + html { + margin: 0; + font: .9em/1.6em "Droid Serif", Cambria, Georgia, "DejaVu Serif", serif; + background-image: url(/img/org-mode-unicorn-logo-worg.png); + background-attachment: fixed; + background-position: right bottom; + background-repeat: no-repeat; + background-color: white; + } + + body { + font-size: 14pt; + line-height: 22pt; + color: black; + margin-top: 0; + } + + body #content { + padding-top: 2em; + margin: auto; + max-width: 70%; + background-color: white; + } + + body #support { + position: fixed; + top:0; + display:block; + font-size: 12pt; + right:0pt; + text-align: right; + padding: .2em 1em; + background: #EEE; + border-radius: 10px; + } + + body .title { + margin-left: 0px; + font-size: 22pt; + } + + #org-div-home-and-up{ + position: fixed; + right: 0.5em; + margin-top: 70px; + font-family:sans-serif; + } + + /* TOC inspired by http://jashkenas.github.com/coffee-script */ + #table-of-contents { + margin-top: 105px; + font-size: 10pt; + font-family:sans-serif; + position: fixed; + right: 0em; + top: 0em; + background: white; + line-height: 12pt; + text-align: right; + box-shadow: 0 0 1em #777777; + -webkit-box-shadow: 0 0 1em #777777; + -moz-box-shadow: 0 0 1em #777777; + -webkit-border-bottom-left-radius: 5px; + -moz-border-radius-bottomleft: 5px; + /* ensure doesn't flow off the screen when expanded */ + max-height: 80%; + overflow: auto; } + #table-of-contents h2 { + font-size: 13pt; + max-width: 9em; + border: 0; + font-weight: normal; + padding-left: 0.5em; + padding-right: 0.5em; + padding-top: 0.05em; + padding-bottom: 0.05em; } + #table-of-contents #text-table-of-contents { + display: none; + text-align: left; } + #table-of-contents:hover #text-table-of-contents { + display: block; + padding: 0.5em; + margin-top: -1.5em; } + + #license { + background-color: #eeeeee; + } + + h1 { + font-size:2.1em; + padding:0 0 30px 0; + margin-top: 10px; + margin-bottom: 10px; + margin-right: 7%; + color: grey; + } + + h2 { + font-family:sans-serif; + font-size:1.45em; + padding:10px 0 10px 0; + color: black; + border-bottom: 1px solid #ddd; + padding-top: 1.5em; + } + + .outline-text-2 { + margin-left: 0.1em + } + + h3 { + font-family:sans-serif; + font-size:1.3em; + color: grey; + margin-left: 0.6em; + padding-top: 1.5em; + } + + /* #A34D32;*/ + + + .outline-text-3 { + margin-left: 0.9em; + } + + h4 { + font-family:sans-serif; + font-size:1.2em; + margin-left: 1.2em; + color: #A5573E; + padding-top: 1.5em; + } + + .outline-text-4 { + margin-left: 1.45em; + } + + a {text-decoration: none; font-weight: 400;} + a:visited {text-decoration: none; font-weight: 400;} + a:hover {text-decoration: underline;} + + .todo { + color: #CA0000; + } + + .done { + color: #006666; + } + + .timestamp-kwd { + color: #444; + } + + .tag { + background-color: #ffff; + color: #ffff; + } + + li { + margin: .4em; + } + + table { + border: 1; + border-color: grey; + } + + thead { + border: 0; + } + + tbody { + border: 0; + } + + tr { + border: 0; + } + + td { + border-left: 0px; + border-right: 0px; + border-top: 0px; + border-bottom: 0px; + } + + th { + border-left: 0px; + border-right: 0px; + border-top: 1px solid grey; + border-bottom: 1px solid grey; + } + + code { + font-size: 100%; + color: black; + padding: 0px 0.2em; + } + + img { + border: 0; + } + + .share img { + opacity: .4; + -moz-opacity: .4; + filter: alpha(opacity=40); + } + + .share img:hover { + opacity: 1; + -moz-opacity: 1; + filter: alpha(opacity=100); + } + + pre { + font-family: Droid Sans Mono, Monaco, Consolas, "Lucida Console", monospace; + color: black; + font-size: 90%; + padding: 0.5em; + overflow: auto; + border: none; + background-color: #f2f2f2; + border-radius: 5px; + } + + .org-info-box { + clear:both; + margin-left:auto; + margin-right:auto; + padding:0.7em; + } + .org-info-box img { + float:left; + margin:0em 0.5em 0em 0em; + } + .org-info-box p { + margin:0em; + padding:0em; + } + + + .builtin { + /* font-lock-builtin-face */ + color: #f4a460; + } + .comment { + /* font-lock-comment-face */ + color: #737373; + } + .comment-delimiter { + /* font-lock-comment-delimiter-face */ + color: #666666; + } + .constant { + /* font-lock-constant-face */ + color: #db7093; + } + .doc { + /* font-lock-doc-face */ + color: #b3b3b3; + } + .function-name { + /* font-lock-function-name-face */ + color: #5f9ea0; + } + .headline { + /* headline-face */ + color: #ffffff; + background-color: #000000; + font-weight: bold; + } + .keyword { + /* font-lock-keyword-face */ + color: #4682b4; + } + .negation-char { + } + .regexp-grouping-backslash { + } + .regexp-grouping-construct { + } + .string { + /* font-lock-string-face */ + color: #ccc79a; + } + .todo-comment { + /* todo-comment-face */ + color: #ffffff; + background-color: #000000; + font-weight: bold; + } + .variable-name { + /* font-lock-variable-name-face */ + color: #ff6a6a; + } + .warning { + /* font-lock-warning-face */ + color: #ffffff; + background-color: #cd5c5c; + font-weight: bold; + } + .important { + /* font-lock-warning-face */ + background-color: #e3e3f7; + } + .exercise { + /* font-lock-warning-face */ + background-color: #e3f7e3; + } + .note { + /* font-lock-warning-face */ + background-color: #f7f7d9; + } + pre.a { + color: inherit; + background-color: inherit; + font: inherit; + text-decoration: inherit; + } + pre.a:hover { + text-decoration: underline; + } + + /* Styles for org-info.js */ + + .org-info-js_info-navigation + { + border-style:none; + } + + #org-info-js_console-label + { + font-size:10px; + font-weight:bold; + white-space:nowrap; + } + + .org-info-js_search-highlight + { + background-color:#ffff00; + color:#000000; + font-weight:bold; + } + + #org-info-js-window + { + border-bottom:1px solid black; + padding-bottom:10px; + margin-bottom:10px; + } + + + + .org-info-search-highlight + { + background-color:#adefef; /* same color as emacs default */ + color:#000000; + font-weight:bold; + } + + .org-bbdb-company { + /* bbdb-company */ + font-style: italic; + } + .org-bbdb-field-name { + } + .org-bbdb-field-value { + } + .org-bbdb-name { + /* bbdb-name */ + text-decoration: underline; + } + .org-bold { + /* bold */ + font-weight: bold; + } + .org-bold-italic { + /* bold-italic */ + font-weight: bold; + font-style: italic; + } + .org-border { + /* border */ + background-color: #000000; + } + .org-buffer-menu-buffer { + /* buffer-menu-buffer */ + font-weight: bold; + } + .org-builtin { + /* font-lock-builtin-face */ + color: #da70d6; + } + .org-button { + /* button */ + text-decoration: underline; + } + .org-c-nonbreakable-space { + /* c-nonbreakable-space-face */ + background-color: #ff0000; + font-weight: bold; + } + .org-calendar-today { + /* calendar-today */ + text-decoration: underline; + } + .org-comment { + /* font-lock-comment-face */ + color: #b22222; + } + .org-comment-delimiter { + /* font-lock-comment-delimiter-face */ + color: #b22222; + } + .org-constant { + /* font-lock-constant-face */ + color: #5f9ea0; + } + .org-cursor { + /* cursor */ + background-color: #000000; + } + .org-default { + /* default */ + color: #000000; + background-color: #ffffff; + } + .org-diary { + /* diary */ + color: #ff0000; + } + .org-doc { + /* font-lock-doc-face */ + color: #bc8f8f; + } + .org-escape-glyph { + /* escape-glyph */ + color: #a52a2a; + } + .org-file-name-shadow { + /* file-name-shadow */ + color: #7f7f7f; + } + .org-fixed-pitch { + } + .org-fringe { + /* fringe */ + background-color: #f2f2f2; + } + .org-function-name { + /* font-lock-function-name-face */ + color: #0000ff; + } + .org-header-line { + /* header-line */ + color: #333333; + background-color: #e5e5e5; + } + .org-help-argument-name { + /* help-argument-name */ + font-style: italic; + } + .org-highlight { + /* highlight */ + background-color: #b4eeb4; + } + .org-holiday { + /* holiday */ + background-color: #ffc0cb; + } + .org-info-header-node { + /* info-header-node */ + color: #a52a2a; + font-weight: bold; + font-style: italic; + } + .org-info-header-xref { + /* info-header-xref */ + color: #0000ff; + text-decoration: underline; + } + .org-info-menu-header { + /* info-menu-header */ + font-weight: bold; + } + .org-info-menu-star { + /* info-menu-star */ + color: #ff0000; + } + .org-info-node { + /* info-node */ + color: #a52a2a; + font-weight: bold; + font-style: italic; + } + .org-info-title-1 { + /* info-title-1 */ + font-size: 172%; + font-weight: bold; + } + .org-info-title-2 { + /* info-title-2 */ + font-size: 144%; + font-weight: bold; + } + .org-info-title-3 { + /* info-title-3 */ + font-size: 120%; + font-weight: bold; + } + .org-info-title-4 { + /* info-title-4 */ + font-weight: bold; + } + .org-info-xref { + /* info-xref */ + color: #0000ff; + text-decoration: underline; + } + .org-isearch { + /* isearch */ + color: #b0e2ff; + background-color: #cd00cd; + } + .org-italic { + /* italic */ + font-style: italic; + } + .org-keyword { + /* font-lock-keyword-face */ + color: #a020f0; + } + .org-lazy-highlight { + /* lazy-highlight */ + background-color: #afeeee; + } + .org-link { + /* link */ + color: #0000ff; + text-decoration: underline; + } + .org-link-visited { + /* link-visited */ + color: #8b008b; + text-decoration: underline; + } + .org-match { + /* match */ + background-color: #ffff00; + } + .org-menu { + } + .org-message-cited-text { + /* message-cited-text */ + color: #ff0000; + } + .org-message-header-cc { + /* message-header-cc */ + color: #191970; + } + .org-message-header-name { + /* message-header-name */ + color: #6495ed; + } + .org-message-header-newsgroups { + /* message-header-newsgroups */ + color: #00008b; + font-weight: bold; + font-style: italic; + } + .org-message-header-other { + /* message-header-other */ + color: #4682b4; + } + .org-message-header-subject { + /* message-header-subject */ + color: #000080; + font-weight: bold; + } + .org-message-header-to { + /* message-header-to */ + color: #191970; + font-weight: bold; + } + .org-message-header-xheader { + /* message-header-xheader */ + color: #0000ff; + } + .org-message-mml { + /* message-mml */ + color: #228b22; + } + .org-message-separator { + /* message-separator */ + color: #a52a2a; + } + .org-minibuffer-prompt { + /* minibuffer-prompt */ + color: #0000cd; + } + .org-mm-uu-extract { + /* mm-uu-extract */ + color: #006400; + background-color: #ffffe0; + } + .org-mode-line { + /* mode-line */ + color: #000000; + background-color: #bfbfbf; + } + .org-mode-line-buffer-id { + /* mode-line-buffer-id */ + font-weight: bold; + } + .org-mode-line-highlight { + } + .org-mode-line-inactive { + /* mode-line-inactive */ + color: #333333; + background-color: #e5e5e5; + } + .org-mouse { + /* mouse */ + background-color: #000000; + } + .org-negation-char { + } + .org-next-error { + /* next-error */ + background-color: #eedc82; + } + .org-nobreak-space { + /* nobreak-space */ + color: #a52a2a; + text-decoration: underline; + } + .org-org-agenda-date { + /* org-agenda-date */ + color: #0000ff; + } + .org-org-agenda-date-weekend { + /* org-agenda-date-weekend */ + color: #0000ff; + font-weight: bold; + } + .org-org-agenda-restriction-lock { + /* org-agenda-restriction-lock */ + background-color: #ffff00; + } + .org-org-agenda-structure { + /* org-agenda-structure */ + color: #0000ff; + } + .org-org-archived { + /* org-archived */ + color: #7f7f7f; + } + .org-org-code { + /* org-code */ + color: #7f7f7f; + } + .org-org-column { + /* org-column */ + background-color: #e5e5e5; + } + .org-org-column-title { + /* org-column-title */ + background-color: #e5e5e5; + font-weight: bold; + text-decoration: underline; + } + .org-org-date { + /* org-date */ + color: #a020f0; + text-decoration: underline; + } + .org-org-done { + /* org-done */ + color: #228b22; + font-weight: bold; + } + .org-org-drawer { + /* org-drawer */ + color: #0000ff; + } + .org-org-ellipsis { + /* org-ellipsis */ + color: #b8860b; + text-decoration: underline; + } + .org-org-formula { + /* org-formula */ + color: #b22222; + } + .org-org-headline-done { + /* org-headline-done */ + color: #bc8f8f; + } + .org-org-hide { + /* org-hide */ + color: #e5e5e5; + } + .org-org-latex-and-export-specials { + /* org-latex-and-export-specials */ + color: #8b4513; + } + .org-org-level-1 { + /* org-level-1 */ + color: #0000ff; + } + .org-org-level-2 { + /* org-level-2 */ + color: #b8860b; + } + .org-org-level-3 { + /* org-level-3 */ + color: #a020f0; + } + .org-org-level-4 { + /* org-level-4 */ + color: #b22222; + } + .org-org-level-5 { + /* org-level-5 */ + color: #228b22; + } + .org-org-level-6 { + /* org-level-6 */ + color: #5f9ea0; + } + .org-org-level-7 { + /* org-level-7 */ + color: #da70d6; + } + .org-org-level-8 { + /* org-level-8 */ + color: #bc8f8f; + } + .org-org-link { + /* org-link */ + color: #a020f0; + text-decoration: underline; + } + .org-org-property-value { + } + .org-org-scheduled-previously { + /* org-scheduled-previously */ + color: #b22222; + } + .org-org-scheduled-today { + /* org-scheduled-today */ + color: #006400; + } + .org-org-sexp-date { + /* org-sexp-date */ + color: #a020f0; + } + .org-org-special-keyword { + /* org-special-keyword */ + color: #bc8f8f; + } + .org-org-table { + /* org-table */ + color: #0000ff; + } + .org-org-tag { + /* org-tag */ + font-weight: bold; + } + .org-org-target { + /* org-target */ + text-decoration: underline; + } + .org-org-time-grid { + /* org-time-grid */ + color: #b8860b; + } + .org-org-todo { + /* org-todo */ + color: #ff0000; + } + .org-org-upcoming-deadline { + /* org-upcoming-deadline */ + color: #b22222; + } + .org-org-verbatim { + /* org-verbatim */ + color: #7f7f7f; + text-decoration: underline; + } + .org-org-warning { + /* org-warning */ + color: #ff0000; + font-weight: bold; + } + .org-outline-1 { + /* outline-1 */ + color: #0000ff; + } + .org-outline-2 { + /* outline-2 */ + color: #b8860b; + } + .org-outline-3 { + /* outline-3 */ + color: #a020f0; + } + .org-outline-4 { + /* outline-4 */ + color: #b22222; + } + .org-outline-5 { + /* outline-5 */ + color: #228b22; + } + .org-outline-6 { + /* outline-6 */ + color: #5f9ea0; + } + .org-outline-7 { + /* outline-7 */ + color: #da70d6; + } + .org-outline-8 { + /* outline-8 */ + color: #bc8f8f; + } + .outline-text-1, .outline-text-2, .outline-text-3, .outline-text-4, .outline-text-5, .outline-text-6 { + /* Add more spacing between section. Padding, so that folding with org-info.js works as expected. */ + + } + + .org-preprocessor { + /* font-lock-preprocessor-face */ + color: #da70d6; + } + .org-query-replace { + /* query-replace */ + color: #b0e2ff; + background-color: #cd00cd; + } + .org-regexp-grouping-backslash { + /* font-lock-regexp-grouping-backslash */ + font-weight: bold; + } + .org-regexp-grouping-construct { + /* font-lock-regexp-grouping-construct */ + font-weight: bold; + } + .org-region { + /* region */ + background-color: #eedc82; + } + .org-rmail-highlight { + } + .org-scroll-bar { + /* scroll-bar */ + background-color: #bfbfbf; + } + .org-secondary-selection { + /* secondary-selection */ + background-color: #ffff00; + } + .org-shadow { + /* shadow */ + color: #7f7f7f; + } + .org-show-paren-match { + /* show-paren-match */ + background-color: #40e0d0; + } + .org-show-paren-mismatch { + /* show-paren-mismatch */ + color: #ffffff; + background-color: #a020f0; + } + .org-string { + /* font-lock-string-face */ + color: #bc8f8f; + } + .org-texinfo-heading { + /* texinfo-heading */ + color: #0000ff; + } + .org-tool-bar { + /* tool-bar */ + color: #000000; + background-color: #bfbfbf; + } + .org-tooltip { + /* tooltip */ + color: #000000; + background-color: #ffffe0; + } + .org-trailing-whitespace { + /* trailing-whitespace */ + background-color: #ff0000; + } + .org-type { + /* font-lock-type-face */ + color: #228b22; + } + .org-underline { + /* underline */ + text-decoration: underline; + } + .org-variable-name { + /* font-lock-variable-name-face */ + color: #b8860b; + } + .org-variable-pitch { + } + .org-vertical-border { + } + .org-warning { + /* font-lock-warning-face */ + color: #ff0000; + font-weight: bold; + } + .rss_box {} + .rss_title, rss_title a {} + .rss_items {} + .rss_item a:link, .rss_item a:visited, .rss_item a:active {} + .rss_item a:hover {} + .rss_date {} + + pre.src { + position: static; + overflow: visible; + padding-top: 1.2em; + } + + label.org-src-name { + font-size: 80%; + font-style: italic; + } + + #show_source {margin: 0; padding: 0;} + + #postamble { + font-size: 75%; + min-width: 700px; + max-width: 80%; + line-height: 14pt; + margin-left: 20px; + margin-top: 10px; + padding: .2em; + background-color: #ffffff; + z-index: -1000; + } + + +} /* END OF @media all */ + +@media screen +{ + #table-of-contents { + position: fixed; + margin-top: 105px; + float: right; + border: 1px solid #red; + max-width: 50%; + overflow: auto; + } +} /* END OF @media screen */ diff --git a/qmckl.html b/qmckl.html new file mode 100644 index 0000000..7673484 --- /dev/null +++ b/qmckl.html @@ -0,0 +1,720 @@ + + + + + + + +Introduction + + + + + + + + + + + + + +
+ UP + | + HOME +
+

Introduction

+ + +
+

1 Using QMCkl

+
+

+The qmckl.h header file has to be included in C codes when +QMCkl functions are used: +

+ +
+
#include "qmckl.h"
+
+
+ +

+In Fortran programs, the qmckl_f.f90 interface file should be +included in the source code using the library, and the Fortran codes +should use the qmckl module as +

+ +
+
use qmckl
+
+
+ +

+Both files are located in the include/ directory. +

+
+
+ +
+

2 Developing in QMCkl

+
+
+
+

2.1 Literate programming

+
+

+In a traditional source code, most of the lines of source files of a program +are code, scripts, Makefiles, and only a few lines are comments explaining +parts of the code that are non-trivial to understand. The documentation of +the prorgam is usually written in a separate directory, and is often outdated +compared to the code. +

+ +

+Literate programming is a different approach to programming, +where the program is considered as a publishable-quality document. Most of +the lines of the source files are text, mathematical formulas, tables, +figures, etc, and the lines of code are just the translation in a computer +language of the ideas and algorithms expressed in the text. More importantly, +the "document" is structured like a text document with sections, subsections, +a bibliography, a table of contents etc, and the place where pieces of code +appear are the places where they should belong for the reader to understand +the logic of the program, not the places where the compiler expects to find +them. Both the publishable-quality document and the binary executable are +produced from the same source files. +

+ +

+Literate programming is particularly well adapted in this context, as the +central part of this project is the documentation of an API. The +implementation of the algorithms is just an expression of the algorithms in a +language that can be compiled, so that the correctness of the algorithms can +be tested. +

+ +

+We have chosen to write the source files in org-mode format, +as any text editor can be used to edit org-mode files. To +produce the documentation, there exists multiple possibilities to convert +org-mode files into different formats such as HTML or PDF. The source code is +easily extracted from the org-mode files invoking the Emacs text editor from +the command-line in the Makefile, and then the produced files are compiled. +Moreover, within the Emacs text editor the source code blocks can be executed +interactively, in the same spirit as Jupyter notebooks. +

+
+
+ + +
+

2.2 Source code editing

+
+

+For a tutorial on literate programming with org-mode, follow this link. +

+ +

+Any text editor can be used to edit org-mode files. For a better +user experience Emacs is recommended. For users hating Emacs, it +is good to know that Emacs can behave like Vim when switched into +``Evil'' mode. +

+ +

+In the tools/init.el file, we provide a minimal Emacs configuration +file for vim users. This file should be copied into .emacs.d/init.el. +

+ +

+For users with a preference for Jupyter notebooks, we also provide the +tools/nb_to_org.sh script can convert jupyter notebooks into org-mode +files. +

+ +

+Note that pandoc can be used to convert multiple markdown formats into +org-mode. +

+
+
+ + +
+

2.3 Choice of the programming language

+
+

+Most of the codes of the TREX CoE are written in Fortran with some scripts in +Bash and Python. Outside of the CoE, Fortran is also important (Casino, Amolqc), +and other important languages used by the community are C and C++ (QMCPack, +QWalk), and Julia is gaining in popularity. The library we design should be +compatible with all of these languages. The QMCkl API has to be compatible +with the C language since libraries with a C-compatible API can be used in +every other language. +

+ +

+High-performance versions of the QMCkl, with the same API, will be rewritten by +the experts in HPC. These optimized libraries will be tuned for specific +architectures, among which we can cite x86 based processors, and GPU +accelerators. Nowadays, the most efficient software tools to take advantage of +low-level features of the processor (intrinsics) and of GPUs are for C++ +developers. It is highly probable that the optimized implementations will be +written in C++, and this is agreement with our choice to make the API +C-compatible. +

+ +

+Fortran is one of the most common languages used by the community, and is simple +enough to make the algorithms readable both by experts in QMC, and experts in +HPC. Hence we propose in this pedagogical implementation of QMCkl to use Fortran +to express the QMC algorithms. As the main languages of the library is C, this +implies that the exposed C functions call the Fortran routine. However, for +internal functions related to system programming, the C language is more natural +than Fortran. +

+ +

+The Fortran source files should provide a C interface using the +iso_c_binding module. The name of the Fortran source files should end with +_f.f90 to be properly handled by the Makefile. The names of the functions +defined in Fortran should be the same as those exposed in the API suffixed by +_f. Fortran interfaces should also be written in the qmckl_f.f90 file. +

+ +

+For more guidelines on using Fortran to generate a C interface, see +this link. +

+
+
+ +
+

2.4 Coding rules

+
+

+The authors should follow the recommendations of the +SEI+CERT+C+Coding+Standard. +

+ +
    +
  • Store a new value in pointers immediately after the memory is +freed
    • +
  • Free dynamically allocated memory when no longer needed
    • +
+
+
+ + +
+

2.5 Design of the library

+
+

+The proposed API should allow the library to: deal with memory transfers +between CPU and accelerators, and to use different levels of floating-point +precision. We chose a multi-layered design with low-level and high-level +functions (see below). +

+
+
+ +
+

2.6 Naming conventions

+
+

+To avoid namespace collisions, we use qmckl_ as a prefix for all exported +functions and variables. All exported header files should have a file name +prefixed with qmckl_. +

+ +

+If the name of the org-mode file is xxx.org, the name of the +produced C files should be xxx.c and xxx.h and the name of the +produced Fortran file should be xxx.f90. +

+ +

+Arrays are in uppercase and scalars are in lowercase. +

+ +

+In the names of the variables and functions, only the singular +form is allowed. +

+
+
+ +
+

2.7 Application programming interface

+
+

+In the C language, the number of bits used by the integer types can change +from one architecture to another one. To circumvent this problem, we choose to +use the integer types defined in <stdint.h> where the number of bits used for +the integers are fixed. +

+ +

+To ensure that the library will be easily usable in any other language +than C, we restrict the data types in the interfaces to the following: +

+
    +
  • 32-bit and 64-bit integers, scalars and and arrays (int32_t and int64_t)
    • +
  • 32-bit and 64-bit floats, scalars and and arrays (float and double)
    • +
  • Pointers are always casted into 64-bit integers, even on legacy 32-bit architectures
    • +
  • ASCII strings are represented as a pointers to character arrays +and terminated by a '\0' character (C convention).
    • +
  • Complex numbers can be represented by an array of 2 floats.
    • +
  • Boolean variables are stored as integers, 1 for true and 0 for false
    • +
  • Floating point variables should be by default
    • +
  • double unless explicitly mentioned
    • +
  • integers used for counting should always be int64_t
    • +
+ +

+To facilitate the use in other languages than C, we will provide some +bindings in other languages in other repositories. +

+
+
+ +
+

2.8 Global state

+
+

+Global variables should be avoided in the library, because it is +possible that one single program needs to use multiple instances +of the library. To solve this problem we propose to use a pointer +to a context variable, built by the library with the +qmckl_context_create function. The =context= contains the global +state of the library, and is used as the first argument of many +QMCkl functions. +

+ +

+The internal structure of the context is not specified, to give a +maximum of freedom to the different implementations. Modifying +the state is done by setters and getters, prefixed by +qmckl_context_set_ an qmckl_context_get_. When a context +variable is modified by a setter, a copy of the old data structure +is made and updated, and the pointer to the new data structure is +returned, such that the old contexts can still be accessed. It is +also possible to modify the state in an impure fashion, using the +qmckl_context_update_ functions. The context and its old +versions can be destroyed with qmckl_context_destroy. +

+
+
+ +
+

2.9 Low-level functions

+
+

+Low-level functions are very simple functions which are leaves of +the function call tree (they don't call any other QMCkl function). +

+ +

+These functions are pure, and unaware of the QMCkl +context. They are not allowed to allocate/deallocate memory, and +if they need temporary memory it should be provided in input. +

+
+
+ +
+

2.10 High-level functions

+
+

+High-level functions are at the top of the function call tree. +They are able to choose which lower-level function to call +depending on the required precision, and do the corresponding type +conversions. These functions are also responsible for allocating +temporary storage, to simplify the use of accelerators. +

+ +

+The high-level functions should be pure, unless the introduction +of non-purity is justified. All the side effects should be made in +the context variable. +

+
+
+ +
+

2.11 Numerical precision

+
+

+The number of bits of precision required for a function should be +given as an input of low-level computational functions. This input +will be used to define the values of the different thresholds that +might be used to avoid computing unnecessary noise. High-level +functions will use the precision specified in the context +variable. +

+
+
+ +
+

2.12 Algorithms

+
+

+Reducing the scaling of an algorithm usually implies also reducing +its arithmetic complexity (number of flops per byte). Therefore, +for small sizes \(\mathcal{O}(N^3)\) and \(\mathcal{O}(N^2)\) +algorithms are better adapted than linear scaling algorithms. As +QMCkl is a general purpose library, multiple algorithms should be +implemented adapted to different problem sizes. +

+
+
+ +
+

2.13 Rules for the API

+
+
    +
  • stdint should be used for integers (int32_t, int64_t)
    • +
  • integers used for counting should always be int64_t
    • +
  • floats should be by default double, unless explicitly mentioned
    • +
  • pointers are converted to int64_t to increase portability
    • +
+
+
+
+
+
+

Author: TREX CoE

+

Created: 2021-03-19 Fri 12:49

+

Validate

+
+ + diff --git a/qmckl_ao.html b/qmckl_ao.html new file mode 100644 index 0000000..3da4bb2 --- /dev/null +++ b/qmckl_ao.html @@ -0,0 +1,1206 @@ + + + + + + + +Atomic Orbitals + + + + + + + + + + + + + +
+ UP + | + HOME +
+

Atomic Orbitals

+ +

+The atomic basis set is defined as a list of shells. Each shell \(s\) is +centered on a nucleus \(A\), possesses a given angular momentum \(l\) and a +radial function \(R_s\). The radial function is a linear combination of +\emph{primitive} functions that can be of type Slater (\(p=1\)) or +Gaussian (\(p=2\)): +

+ +

+\[ + R_s(\mathbf{r}) = \mathcal{N}_s |\mathbf{r}-\mathbf{R}_A|^{n_s} + \sum_{k=1}^{N_{\text{prim}}} a_{ks} + \exp \left( - \gamma_{ks} | \mathbf{r}-\mathbf{R}_A | ^p \right). | +\] +

+ +

+In the case of Gaussian functions, \(n_s\) is always zero. +The normalization factor \(\mathcal{N}_s\) ensures that all the functions +of the shell are normalized to unity. As this normalization requires +the ability to compute overlap integrals, it should be written in the +file to ensure that the file is self-contained and does not require +the client program to have the ability to compute such integrals. +

+ +

+Atomic orbitals (AOs) are defined as +

+ +

+\[ +\chi_i (\mathbf{r}) = P_{\eta(i)}(\mathbf{r})\, R_{\theta(i)} (\mathbf{r}) +\] +

+ +

+where \(\theta(i)\) returns the shell on which the AO is expanded, +and \(\eta(i)\) denotes which angular function is chosen. +

+ +

+In this section we describe the kernels used to compute the values, +gradients and Laplacian of the atomic basis functions. +

+ +
+

1 Polynomial part

+
+
+
+

1.1 Powers of \(x-X_i\)

+
+

+The qmckl_ao_power function computes all the powers of the n +input data up to the given maximum value given in input for each of +the \(n\) points: +

+ +

+\[ P_{ik} = X_i^k \] +

+ + + + +++ ++ ++ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
contextinputGlobal state
ninputNumber of values
X(n)inputArray containing the input values
LMAX(n)inputArray containing the maximum power for each value
P(LDP,n)outputArray containing all the powers of X
LDPinputLeading dimension of array P
+ +

+Requirements: +

+ +
    +
  • context is not QMCKL_NULL_CONTEXT
    • +
  • n > 0
    • +
  • X is allocated with at least \(n \times 8\) bytes
    • +
  • LMAX is allocated with at least \(n \times 4\) bytes
    • +
  • P is allocated with at least \(n \times \max_i \text{LMAX}_i \times 8\) bytes
    • +
  • LDP >= \(\max_i\) LMAX[i]
    • +
+ +
+
qmckl_exit_code
+qmckl_ao_power(const qmckl_context context,
+               const int64_t n, 
+               const double *X,
+               const int32_t *LMAX,
+               const double *P,
+               const int64_t LDP);
+
+
+ +
+
integer function qmckl_ao_power_f(context, n, X, LMAX, P, ldp) result(info)
+  use qmckl
+  implicit none
+  integer*8 , intent(in)  :: context
+  integer*8 , intent(in)  :: n
+  real*8    , intent(in)  :: X(n)
+  integer   , intent(in)  :: LMAX(n)
+  real*8    , intent(out) :: P(ldp,n)
+  integer*8 , intent(in)  :: ldp
+
+  integer*8  :: i,k
+
+  info = QMCKL_SUCCESS
+
+  if (context == QMCKL_NULL_CONTEXT) then
+     info = QMCKL_INVALID_CONTEXT
+     return
+  endif
+
+  if (n <= ldp) then
+     info = QMCKL_INVALID_ARG_2
+     return
+  endif
+
+  k = MAXVAL(LMAX)
+  if (LDP < k) then
+     info = QMCKL_INVALID_ARG_6
+     return
+  endif
+
+  if (k <= 0) then
+     info = QMCKL_INVALID_ARG_4
+     return
+  endif
+
+  do i=1,n
+     P(1,i) = X(i)
+     do k=2,LMAX(i)
+        P(k,i) = P(k-1,i) * X(i) 
+     end do
+  end do
+
+end function qmckl_ao_power_f
+
+
+ +
+
integer(c_int32_t) function test_qmckl_ao_power(context) bind(C)
+  use qmckl
+  implicit none
+
+  integer(c_int64_t), intent(in), value :: context
+
+  integer*8                     :: n, LDP 
+  integer, allocatable          :: LMAX(:) 
+  double precision, allocatable :: X(:), P(:,:)
+  integer*8                     :: i,j
+  double precision              :: epsilon
+
+  epsilon = qmckl_context_get_epsilon(context)
+
+  n = 100;
+  LDP = 10;
+
+  allocate(X(n), P(LDP,n), LMAX(n))
+
+  do j=1,n
+     X(j) = -5.d0 + 0.1d0 * dble(j)
+     LMAX(j) = 1 + int(mod(j, 5),4)
+  end do
+
+  test_qmckl_ao_power = qmckl_ao_power(context, n, X, LMAX, P, LDP) 
+  if (test_qmckl_ao_power /= 0) return
+
+  test_qmckl_ao_power = -1
+
+  do j=1,n
+     do i=1,LMAX(j)
+        if ( X(j)**i == 0.d0 ) then
+           if ( P(i,j) /= 0.d0) return
+        else
+           if ( dabs(1.d0 - P(i,j) / (X(j)**i)) > epsilon ) return
+        end if
+     end do
+  end do
+
+  test_qmckl_ao_power = 0
+  deallocate(X,P,LMAX)
+end function test_qmckl_ao_power
+
+
+
+
+ +
+

1.2 Value, Gradient and Laplacian of a polynomial

+
+

+A polynomial is centered on a nucleus \(\mathbf{R}_i\) +

+ +

+\[ + P_l(\mathbf{r},\mathbf{R}_i) = (x-X_i)^a (y-Y_i)^b (z-Z_i)^c + \] +

+ +

+The gradients with respect to electron coordinates are +

+ +\begin{eqnarray*} +\frac{\partial }{\partial x} P_l\left(\mathbf{r},\mathbf{R}_i \right) & + = & a (x-X_i)^{a-1} (y-Y_i)^b (z-Z_i)^c \\ +\frac{\partial }{\partial y} P_l\left(\mathbf{r},\mathbf{R}_i \right) & + = & b (x-X_i)^a (y-Y_i)^{b-1} (z-Z_i)^c \\ +\frac{\partial }{\partial z} P_l\left(\mathbf{r},\mathbf{R}_i \right) & + = & c (x-X_i)^a (y-Y_i)^b (z-Z_i)^{c-1} \\ +\end{eqnarray*} + +

+and the Laplacian is +

+ +\begin{eqnarray*} +\left( \frac{\partial }{\partial x^2} + + \frac{\partial }{\partial y^2} + + \frac{\partial }{\partial z^2} \right) P_l + \left(\mathbf{r},\mathbf{R}_i \right) & = & + a(a-1) (x-X_i)^{a-2} (y-Y_i)^b (z-Z_i)^c + \\ + && b(b-1) (x-X_i)^a (y-Y_i)^{b-1} (z-Z_i)^c + \\ + && c(c-1) (x-X_i)^a (y-Y_i)^b (z-Z_i)^{c-1}. +\end{eqnarray*} + +

+qmckl_ao_polynomial_vgl computes the values, gradients and +Laplacians at a given point in space, of all polynomials with an +angular momentum up to lmax. +

+ + + + +++ ++ ++ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
contextinputGlobal state
X(3)inputArray containing the coordinates of the points
R(3)inputArray containing the x,y,z coordinates of the center
lmaxinputMaximum angular momentum
noutputNumber of computed polynomials
L(ldl,n)outputContains a,b,c for all n results
ldlinputLeading dimension of L
VGL(ldv,n)outputValue, gradients and Laplacian of the polynomials
ldvinputLeading dimension of array VGL
+ +

+Requirements: +

+ +
    +
  • context is not QMCKL_NULL_CONTEXT
    • +
  • n > 0
    • +
  • lmax >= 0
    • +
  • ldl >= 3
    • +
  • ldv >= 5
    • +
  • X is allocated with at least \(3 \times 8\) bytes
    • +
  • R is allocated with at least \(3 \times 8\) bytes
    • +
  • n >= (lmax+1)(lmax+2)(lmax+3)/6
    • +
  • L is allocated with at least \(3 \times n \times 4\) bytes
    • +
  • VGL is allocated with at least \(5 \times n \times 8\) bytes
    • +
  • On output, n should be equal to (lmax+1)(lmax+2)(lmax+3)/6
    • +
  • On output, the powers are given in the following order (l=a+b+c): +
      +
    • Increasing values of l
      • +
    • Within a given value of l, alphabetical order of the +string made by a*"x" + b*"y" + c*"z" (in Python notation). +For example, with a=0, b=2 and c=1 the string is "yyz"
      • +
    • +
+ +
+
qmckl_exit_code
+qmckl_ao_polynomial_vgl(const qmckl_context context,
+                        const double *X,
+                        const double *R,
+                        const int32_t lmax,
+                        const int64_t *n,
+                        const int32_t *L,
+                        const int64_t ldl,
+                        const double *VGL,
+                        const int64_t ldv);
+
+
+ +
+
integer function qmckl_ao_polynomial_vgl_f(context, X, R, lmax, n, L, ldl, VGL, ldv) result(info)
+  use qmckl
+  implicit none
+  integer*8 , intent(in)  :: context
+  real*8    , intent(in)  :: X(3), R(3)
+  integer   , intent(in)  :: lmax
+  integer*8 , intent(out) :: n
+  integer   , intent(out) :: L(ldl,(lmax+1)*(lmax+2)*(lmax+3)/6)
+  integer*8 , intent(in)  :: ldl
+  real*8    , intent(out) :: VGL(ldv,(lmax+1)*(lmax+2)*(lmax+3)/6)
+  integer*8 , intent(in)  :: ldv
+
+  integer*8         :: i,j
+  integer           :: a,b,c,d
+  real*8            :: Y(3)
+  integer           :: lmax_array(3)
+  real*8            :: pows(-2:lmax,3)
+  integer, external :: qmckl_ao_power_f
+  double precision  :: xy, yz, xz
+  double precision  :: da, db, dc, dd
+
+  info = 0
+
+  if (context == QMCKL_NULL_CONTEXT) then
+     info = QMCKL_INVALID_CONTEXT
+     return
+  endif
+
+  if (lmax <= 0) then
+     info = QMCKL_INVALID_ARG_4
+     return
+  endif
+
+  if (ldl < 3) then
+     info = QMCKL_INVALID_ARG_7
+     return
+  endif
+
+  if (ldv < 5) then
+     info = QMCKL_INVALID_ARG_9
+     return
+  endif
+
+
+  do i=1,3
+     Y(i) = X(i) - R(i)
+  end do
+
+  lmax_array(1:3) = lmax
+  if (lmax == 0) then
+     VGL(1,1) = 1.d0
+     vgL(2:5,1) = 0.d0
+     l(1:3,1) = 0
+     n=1
+  else if (lmax > 0) then
+     pows(-2:0,1:3) = 1.d0
+     do i=1,lmax
+        pows(i,1) = pows(i-1,1) * Y(1) 
+        pows(i,2) = pows(i-1,2) * Y(2) 
+        pows(i,3) = pows(i-1,3) * Y(3) 
+     end do
+
+     VGL(1:5,1:4) = 0.d0
+     l  (1:3,1:4) = 0
+
+     VGL(1  ,1  ) = 1.d0
+     vgl(1:5,2:4) = 0.d0
+
+     l  (1,2) = 1
+     vgl(1,2) = pows(1,1)
+     vgL(2,2) = 1.d0
+
+     l  (2,3) = 1
+     vgl(1,3) = pows(1,2)
+     vgL(3,3) = 1.d0
+
+     l  (3,4) = 1
+     vgl(1,4) = pows(1,3)
+     vgL(4,4) = 1.d0
+
+     n=4
+  endif
+
+  ! l>=2
+  dd = 2.d0
+  do d=2,lmax
+     da = dd
+     do a=d,0,-1
+        db = dd-da
+        do b=d-a,0,-1
+           c  = d  - a  - b
+           dc = dd - da - db
+           n = n+1
+
+           l(1,n) = a
+           l(2,n) = b
+           l(3,n) = c
+
+           xy = pows(a,1) * pows(b,2)
+           yz = pows(b,2) * pows(c,3)
+           xz = pows(a,1) * pows(c,3)
+
+           vgl(1,n) = xy * pows(c,3)
+
+           xy = dc * xy
+           xz = db * xz
+           yz = da * yz
+
+           vgl(2,n) = pows(a-1,1) * yz
+           vgl(3,n) = pows(b-1,2) * xz
+           vgl(4,n) = pows(c-1,3) * xy
+
+           vgl(5,n) = &
+                (da-1.d0) * pows(a-2,1) * yz + &
+                (db-1.d0) * pows(b-2,2) * xz + &
+                (dc-1.d0) * pows(c-2,3) * xy
+
+           db = db - 1.d0
+        end do
+        da = da - 1.d0
+     end do
+     dd = dd + 1.d0
+  end do
+
+  info = QMCKL_SUCCESS
+
+end function qmckl_ao_polynomial_vgl_f
+
+
+ + +
+
integer(c_int32_t) function test_qmckl_ao_polynomial_vgl(context) bind(C)
+  use qmckl
+  implicit none
+
+  integer(c_int64_t), intent(in), value :: context
+
+  integer                       :: lmax, d, i
+  integer, allocatable          :: L(:,:)
+  integer*8                     :: n, ldl, ldv, j
+  double precision              :: X(3), R(3), Y(3)
+  double precision, allocatable :: VGL(:,:)
+  double precision              :: w
+  double precision              :: epsilon
+
+  epsilon = qmckl_context_get_epsilon(context)
+
+  X = (/ 1.1 , 2.2 ,  3.3 /)
+  R = (/ 0.1 , 1.2 , -2.3 /)
+  Y(:) = X(:) - R(:)
+
+  lmax = 4;
+  ldl = 3;
+  ldv = 100;
+
+  d = (lmax+1)*(lmax+2)*(lmax+3)/6
+
+  allocate (L(ldl,d), VGL(ldv,d))
+
+  test_qmckl_ao_polynomial_vgl = &
+       qmckl_ao_polynomial_vgl(context, X, R, lmax, n, L, ldl, VGL, ldv)
+
+  if (test_qmckl_ao_polynomial_vgl /= QMCKL_SUCCESS) return
+  if (n /= d) return
+
+  do j=1,n
+     test_qmckl_ao_polynomial_vgl = QMCKL_FAILURE
+     do i=1,3
+        if (L(i,j) < 0) return
+     end do
+     test_qmckl_ao_polynomial_vgl = QMCKL_FAILURE
+     if (dabs(1.d0 - VGL(1,j) / (&
+          Y(1)**L(1,j) * Y(2)**L(2,j) * Y(3)**L(3,j)  &
+          )) > epsilon ) return
+
+     test_qmckl_ao_polynomial_vgl = QMCKL_FAILURE
+     if (L(1,j) < 1) then
+        if (VGL(2,j) /= 0.d0) return
+     else
+        if (dabs(1.d0 - VGL(2,j) / (&
+             L(1,j) * Y(1)**(L(1,j)-1) * Y(2)**L(2,j) * Y(3)**L(3,j) &
+             )) > epsilon ) return
+     end if
+
+     test_qmckl_ao_polynomial_vgl = QMCKL_FAILURE
+     if (L(2,j) < 1) then
+        if (VGL(3,j) /= 0.d0) return
+     else
+        if (dabs(1.d0 - VGL(3,j) / (&
+             L(2,j) * Y(1)**L(1,j) * Y(2)**(L(2,j)-1) * Y(3)**L(3,j) &
+             )) > epsilon ) return
+     end if
+
+     test_qmckl_ao_polynomial_vgl = QMCKL_FAILURE
+     if (L(3,j) < 1) then
+        if (VGL(4,j) /= 0.d0) return
+     else
+        if (dabs(1.d0 - VGL(4,j) / (&
+             L(3,j) * Y(1)**L(1,j) * Y(2)**L(2,j) * Y(3)**(L(3,j)-1) &
+             )) > epsilon ) return
+     end if
+
+     test_qmckl_ao_polynomial_vgl = QMCKL_FAILURE
+     w = 0.d0
+     if (L(1,j) > 1) then
+        w = w + L(1,j) * (L(1,j)-1) * Y(1)**(L(1,j)-2) * Y(2)**L(2,j) * Y(3)**L(3,j) 
+     end if
+     if (L(2,j) > 1) then
+        w = w + L(2,j) * (L(2,j)-1) * Y(1)**L(1,j) * Y(2)**(L(2,j)-2) * Y(3)**L(3,j) 
+     end if
+     if (L(3,j) > 1) then
+        w = w + L(3,j) * (L(3,j)-1) * Y(1)**L(1,j) * Y(2)**L(2,j) * Y(3)**(L(3,j)-2) 
+     end if
+     if (dabs(1.d0 - VGL(5,j) / w) > epsilon ) return
+  end do
+
+  test_qmckl_ao_polynomial_vgl = QMCKL_SUCCESS
+
+  deallocate(L,VGL)
+end function test_qmckl_ao_polynomial_vgl
+
+
+ +
+
int  test_qmckl_ao_polynomial_vgl(qmckl_context context);
+munit_assert_int(0, ==, test_qmckl_ao_polynomial_vgl(context));
+
+
+
+
+
+ +
+

2 Gaussian basis functions

+
+

+qmckl_ao_gaussian_vgl computes the values, gradients and +Laplacians at a given point of n Gaussian functions centered at +the same point: +

+ +

+\[ v_i = \exp(-a_i |X-R|^2) \] +\[ \nabla_x v_i = -2 a_i (X_x - R_x) v_i \] +\[ \nabla_y v_i = -2 a_i (X_y - R_y) v_i \] +\[ \nabla_z v_i = -2 a_i (X_z - R_z) v_i \] +\[ \Delta v_i = a_i (4 |X-R|^2 a_i - 6) v_i \] +

+ + + + +++ ++ ++ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
contextinputGlobal state
X(3)inputArray containing the coordinates of the points
R(3)inputArray containing the x,y,z coordinates of the center
ninputNumber of computed gaussians
A(n)inputExponents of the Gaussians
VGL(ldv,5)outputValue, gradients and Laplacian of the Gaussians
ldvinputLeading dimension of array VGL
+ +

+Requirements : +

+ +
    +
  • context is not 0
    • +
  • n > 0
    • +
  • ldv >= 5
    • +
  • A(i) > 0 for all i
    • +
  • X is allocated with at least \(3 \times 8\) bytes
    • +
  • R is allocated with at least \(3 \times 8\) bytes
    • +
  • A is allocated with at least \(n \times 8\) bytes
    • +
  • VGL is allocated with at least \(n \times 5 \times 8\) bytes
    • +
+ +
+
qmckl_exit_code
+qmckl_ao_gaussian_vgl(const qmckl_context context,
+                      const double *X,
+                      const double *R,
+                      const int64_t *n,
+                      const int64_t *A,
+                      const double *VGL,
+                      const int64_t ldv);
+
+
+ +
+
integer function qmckl_ao_gaussian_vgl_f(context, X, R, n, A, VGL, ldv) result(info)
+  use qmckl
+  implicit none
+  integer*8 , intent(in)  :: context
+  real*8    , intent(in)  :: X(3), R(3)
+  integer*8 , intent(in)  :: n
+  real*8    , intent(in)  :: A(n)
+  real*8    , intent(out) :: VGL(ldv,5)
+  integer*8 , intent(in)  :: ldv
+
+  integer*8         :: i,j
+  real*8            :: Y(3), r2, t, u, v
+
+  info = QMCKL_SUCCESS
+
+  if (context == QMCKL_NULL_CONTEXT) then
+     info = QMCKL_INVALID_CONTEXT
+     return
+  endif
+
+  if (n <= 0) then
+     info = QMCKL_INVALID_ARG_4
+     return
+  endif
+
+  if (ldv < n) then
+     info = QMCKL_INVALID_ARG_7
+     return
+  endif
+
+
+  do i=1,3
+     Y(i) = X(i) - R(i)
+  end do
+  r2 = Y(1)*Y(1) + Y(2)*Y(2) + Y(3)*Y(3)
+
+  do i=1,n
+     VGL(i,1) = dexp(-A(i) * r2)
+  end do
+
+  do i=1,n
+     VGL(i,5) = A(i) * VGL(i,1)
+  end do
+
+  t = -2.d0 * ( X(1) - R(1) )
+  u = -2.d0 * ( X(2) - R(2) )
+  v = -2.d0 * ( X(3) - R(3) )
+
+  do i=1,n
+     VGL(i,2) = t * VGL(i,5)
+     VGL(i,3) = u * VGL(i,5)
+     VGL(i,4) = v * VGL(i,5)
+  end do
+
+  t = 4.d0 * r2
+  do i=1,n
+     VGL(i,5) = (t * A(i) - 6.d0) *  VGL(i,5)
+  end do
+
+end function qmckl_ao_gaussian_vgl_f
+
+
+ +
+
integer(c_int32_t) function test_qmckl_ao_gaussian_vgl(context) bind(C)
+  use qmckl
+  implicit none
+
+  integer(c_int64_t), intent(in), value :: context
+
+  integer*8                     :: n, ldv, j, i
+  double precision              :: X(3), R(3), Y(3), r2
+  double precision, allocatable :: VGL(:,:), A(:)
+  double precision              :: epsilon
+
+  epsilon = qmckl_context_get_epsilon(context)
+
+  X = (/ 1.1 , 2.2 ,  3.3 /)
+  R = (/ 0.1 , 1.2 , -2.3 /)
+  Y(:) = X(:) - R(:)
+  r2 = Y(1)**2 + Y(2)**2 + Y(3)**2
+
+  n = 10;
+  ldv = 100;
+
+  allocate (A(n), VGL(ldv,5))
+  do i=1,n
+     A(i) = 0.0013 * dble(ishft(1,i))
+  end do
+
+
+  test_qmckl_ao_gaussian_vgl = &
+       qmckl_ao_gaussian_vgl(context, X, R, n, A, VGL, ldv)
+  if (test_qmckl_ao_gaussian_vgl /= 0) return
+
+  test_qmckl_ao_gaussian_vgl = -1
+
+  do i=1,n
+     test_qmckl_ao_gaussian_vgl = -11
+     if (dabs(1.d0 - VGL(i,1) / (&
+          dexp(-A(i) * r2) &
+          )) > epsilon ) return
+
+     test_qmckl_ao_gaussian_vgl = -12
+     if (dabs(1.d0 - VGL(i,2) / (&
+          -2.d0 * A(i) * Y(1) * dexp(-A(i) * r2) &
+          )) > epsilon ) return
+
+     test_qmckl_ao_gaussian_vgl = -13
+     if (dabs(1.d0 - VGL(i,3) / (&
+          -2.d0 * A(i) * Y(2) * dexp(-A(i) * r2) &
+          )) > epsilon ) return
+
+     test_qmckl_ao_gaussian_vgl = -14
+     if (dabs(1.d0 - VGL(i,4) / (&
+          -2.d0 * A(i) * Y(3) * dexp(-A(i) * r2) &
+          )) > epsilon ) return
+
+     test_qmckl_ao_gaussian_vgl = -15
+     if (dabs(1.d0 - VGL(i,5) / (&
+          A(i) * (4.d0*r2*A(i) - 6.d0) * dexp(-A(i) * r2) &
+          )) > epsilon ) return
+  end do
+
+  test_qmckl_ao_gaussian_vgl = 0
+
+  deallocate(VGL)
+end function test_qmckl_ao_gaussian_vgl
+
+
+
+
+ +
+

3 TODO Slater basis functions

+
+
+
+

Author: TREX CoE

+

Created: 2021-03-19 Fri 12:49

+

Validate

+
+ + diff --git a/qmckl_context.html b/qmckl_context.html new file mode 100644 index 0000000..a494ef0 --- /dev/null +++ b/qmckl_context.html @@ -0,0 +1,1869 @@ + + + + + + + +Context + + + + + + + + + + + + + +
+ UP + | + HOME +
+

Context

+ +

+The context variable is a handle for the state of the library, +and is stored in a data structure which can't be seen outside of +the library. To simplify compatibility with other languages, the +pointer to the internal data structure is converted into a 64-bit +signed integer, defined in the qmckl_context type. +A value of QMCKL_NULL_CONTEXT for the context is equivalent to a +NULL pointer. +

+ +
+
typedef int64_t qmckl_context ;
+#define QMCKL_NULL_CONTEXT (qmckl_context) 0
+
+
+ +
+

1 Context handling

+
+

+The context appears as an immutable data structure: modifying a +context returns a new context with the modifications. Therefore, it +is necessary to store a pointer to the old version of context so +that it can be freed when necessary. +Note that we also provide a possibility to mutate the context, but +this should be done with caution, only when it is justified. +

+ +

+By convention, in this file context is a qmckl_context variable +and ctx is a qmckl_context_struct* pointer. +

+
+ +
+

1.1 Data structure

+
+

+The main data structure contains pointers to other data structures, +containing the data specific to each given domain, such that the +modified contexts don't need to duplicate the data but only the +pointers. +

+ +
+
typedef struct qmckl_context_struct {
+
+  /* Pointer to the previous context, before modification */
+  struct qmckl_context_struct * prev;
+
+  /* Molecular system */
+  qmckl_ao_basis_struct    * ao_basis;
+
+  /* To be implemented:
+  qmckl_nucleus_struct     * nucleus;
+  qmckl_electron_struct    * electron;
+  qmckl_mo_struct          * mo;
+  qmckl_determinant_struct * det;
+  */
+
+  /* Numerical precision */
+  qmckl_precision_struct   * fp;
+
+  /* Error handling */
+  qmckl_error_struct       * error;
+
+  /* Memory allocation */
+  qmckl_memory_struct      * alloc;
+
+  /* Thread lock */
+  int                        lock_count;
+  pthread_mutex_t            mutex;
+
+  /* Validity checking */
+  uint32_t                   tag;
+
+} qmckl_context_struct;
+
+
+ + +

+A tag is used internally to check if the memory domain pointed +by a pointer is a valid context. This allows to check that even if +the pointer associated with a context is non-null, we can still +verify that it points to the expected data structure. +

+ +
+
#define VALID_TAG   0xBEEFFACE
+#define INVALID_TAG 0xDEADBEEF
+
+
+ +

+The qmckl_context_check function checks if the domain pointed by +the pointer is a valid context. It returns the input qmckl_context +if the context is valid, QMCKL_NULL_CONTEXT otherwise. +

+ +
+
qmckl_context qmckl_context_check(const qmckl_context context) ;
+
+
+ +
+
qmckl_context qmckl_context_check(const qmckl_context context) {
+
+  if (context == QMCKL_NULL_CONTEXT)
+    return QMCKL_NULL_CONTEXT;
+
+  const qmckl_context_struct* ctx = (qmckl_context_struct*) context;
+
+  /* Try to access memory */
+  if (ctx->tag != VALID_TAG) {
+      return QMCKL_NULL_CONTEXT;
+  }
+
+  return context;
+}
+
+
+
+
+ +
+

1.2 Creation

+
+

+To create a new context, qmckl_context_create() should be used. +

+
    +
  • Upon success, it returns a pointer to a new context with the qmckl_context type
    • +
  • It returns QMCKL_NULL_CONTEXT upon failure to allocate the internal data structure
    • +
+ +
+
qmckl_context qmckl_context_create() {
+
+  qmckl_context_struct* ctx =
+    (qmckl_context_struct*) qmckl_malloc (QMCKL_NULL_CONTEXT, sizeof(qmckl_context_struct));
+
+  if (ctx == NULL) {
+    return QMCKL_NULL_CONTEXT;
+  }
+
+  /* Set all pointers to NULL */
+  memset(ctx, 0, sizeof(qmckl_context_struct));
+
+  /* Initialize lock */
+  init_lock(&(ctx->mutex));
+
+  /* Initialize data */
+  ctx->tag = VALID_TAG;
+
+  const qmckl_context context = (qmckl_context) ctx;
+  assert ( qmckl_context_check(context) != QMCKL_NULL_CONTEXT );
+
+  /*
+  qmckl_memory_struct* alloc = (qmckl_memory_struct*)
+    malloc(sizeof(qmckl_memory_struct));
+
+  if (alloc == NULL) {
+    qmckl_unlock(context);
+    return QMCKL_NULL_CONTEXT;
+  }
+
+  memset(alloc,0,sizeof(qmckl_memory_struct));
+
+  ctx->alloc = alloc;
+  */
+
+  return context;
+}
+
+
+
+
+
+

1.3 Access to the previous context

+
+

+qmckl_context_previous returns the previous context. It returns +QMCKL_NULL_CONTEXT for the initial context and for the NULL context. +

+ +
+
qmckl_context qmckl_context_previous(const qmckl_context context) {
+
+  const qmckl_context checked_context = qmckl_context_check(context);
+  if (checked_context == (qmckl_context) 0) {
+    return (qmckl_context) 0;
+  }
+
+  const qmckl_context_struct* ctx = (qmckl_context_struct*) checked_context;
+  return qmckl_context_check((qmckl_context) ctx->prev);
+}
+
+
+
+
+
+

1.4 Locking

+
+

+For thread safety, the context may be locked/unlocked. The lock is +initialized with the PTHREAD_MUTEX_RECURSIVE attribute, and the +number of times the thread has locked it is saved in the +lock_count attribute. +

+ +
+
void init_lock(pthread_mutex_t* mutex) {
+  pthread_mutexattr_t attr;
+  int rc;
+
+  rc = pthread_mutexattr_init(&attr); 
+  assert (rc == 0);
+
+  (void) pthread_mutexattr_settype(&attr, PTHREAD_MUTEX_RECURSIVE);
+
+  rc = pthread_mutex_init ( mutex, &attr);
+  assert (rc == 0);
+
+  (void)pthread_mutexattr_destroy(&attr);
+}
+
+void qmckl_lock(qmckl_context context) {
+  if (context == QMCKL_NULL_CONTEXT)
+    return ;
+  qmckl_context_struct *ctx = (qmckl_context_struct*) context;
+  errno = 0;
+  int rc = pthread_mutex_lock( &(ctx->mutex) );
+  if (rc != 0) {
+    fprintf(stderr, "qmckl_lock:%s\n", strerror(rc) );
+    fflush(stderr);
+  }
+  assert (rc == 0);
+  ctx->lock_count++;
+/*
+  printf("  lock : %d\n", ctx->lock_count);
+*/
+}
+
+void qmckl_unlock(qmckl_context context) {
+  qmckl_context_struct *ctx = (qmckl_context_struct*) context;
+  int rc = pthread_mutex_unlock( &(ctx->mutex) );
+  if (rc != 0) {
+    fprintf(stderr, "qmckl_unlock:%s\n", strerror(rc) );
+    fflush(stderr);
+  }
+  assert (rc == 0);
+  ctx->lock_count--;
+/*
+  printf("unlock : %d\n", ctx->lock_count);
+*/
+}
+
+
+
+
+ +
+

1.5 Copy

+
+

+qmckl_context_copy makes a shallow copy of a context. It returns +QMCKL_NULL_CONTEXT upon failure. +

+ +
+
qmckl_context qmckl_context_copy(const qmckl_context context) {
+
+  qmckl_lock(context);
+
+  const qmckl_context checked_context = qmckl_context_check(context);
+
+  if (checked_context == QMCKL_NULL_CONTEXT) {
+    qmckl_unlock(context);
+    return QMCKL_NULL_CONTEXT;
+  }
+
+
+  qmckl_context_struct* old_ctx =
+    (qmckl_context_struct*) checked_context;
+
+  qmckl_context_struct* new_ctx =
+    (qmckl_context_struct*) qmckl_malloc (context, sizeof(qmckl_context_struct));
+
+  if (new_ctx == NULL) {
+    qmckl_unlock(context);
+    return QMCKL_NULL_CONTEXT;
+  }
+
+  /* Copy the old context on the new one */
+  memcpy(new_ctx, old_ctx, sizeof(qmckl_context_struct));
+
+  new_ctx->prev = old_ctx;
+
+  qmckl_unlock( (qmckl_context) new_ctx );
+  qmckl_unlock( (qmckl_context) old_ctx );
+
+  return (qmckl_context) new_ctx;
+}
+
+
+
+
+
+
+

1.6 Destroy

+
+

+The context is destroyed with qmckl_context_destroy, leaving the ancestors untouched. +It frees the context, and returns the previous context. +

+ +
+
qmckl_context qmckl_context_destroy(const qmckl_context context) {
+
+  const qmckl_context checked_context = qmckl_context_check(context);
+  if (checked_context == QMCKL_NULL_CONTEXT) return QMCKL_NULL_CONTEXT;
+
+  qmckl_lock(context);
+
+  qmckl_context_struct* ctx = (qmckl_context_struct*) context;
+  assert (ctx != NULL);  /* Shouldn't be true because the context is valid */
+
+  qmckl_unlock(context);
+
+  const qmckl_context prev_context = (qmckl_context) ctx->prev;
+  if (prev_context == QMCKL_NULL_CONTEXT) {
+    /* This is the first context, free all memory. */
+    struct qmckl_memory_struct* new = NULL;
+    while (ctx->alloc != NULL) {
+       new = ctx->alloc->next;
+       free(ctx->alloc->pointer);
+       ctx->alloc->pointer = NULL;
+       free(ctx->alloc);
+       ctx->alloc = new;
+    }
+  }
+
+  qmckl_exit_code rc;
+  rc = qmckl_context_remove_memory(context,ctx);
+/*
+  assert (rc == QMCKL_SUCCESS);
+  */
+
+  ctx->tag = INVALID_TAG;
+
+  const int rc_destroy = pthread_mutex_destroy( &(ctx->mutex) );
+  if (rc_destroy != 0) {
+     fprintf(stderr, "qmckl_context_destroy: %s (count = %d)\n", strerror(rc_destroy), ctx->lock_count);
+     abort();
+  }
+
+  rc = qmckl_free(context,ctx);
+  assert (rc == QMCKL_SUCCESS);
+
+  //memset(ctx, 0, sizeof(qmckl_context_struct));
+
+
+  return prev_context;
+}
+
+
+
+
+
+
+

2 Memory allocation handling

+
+
+
+

2.1 Data structure

+
+

+Pointers to all allocated memory domains are stored in the context, +in a linked list. The size is also stored, to enable the +computation of the amount of currently used memory by the library. +

+ +
+
typedef struct qmckl_memory_struct {
+  struct qmckl_memory_struct * next    ;
+  void                       * pointer ;
+  size_t                       size    ;
+} qmckl_memory_struct;
+
+
+
+
+ +
+

2.2 Append memory

+
+

+The following function, called in qmckl_memory.c, appends a new +pair (pointer, size) to the data structure. +It is forbidden to pass the NULL pointer, or a zero size. +If the context is QMCKL_NULL_CONTEXT, the function returns +immediately with QMCKL_SUCCESS. +

+ +
+
qmckl_exit_code qmckl_context_append_memory(qmckl_context context,
+                                            void* pointer,
+                                            const size_t size) {
+  assert (pointer != NULL);
+  assert (size > 0L);
+
+  qmckl_lock(context);
+
+  if ( qmckl_context_check(context) == QMCKL_NULL_CONTEXT ) {
+    qmckl_unlock(context);
+    return QMCKL_SUCCESS;
+  }
+
+  qmckl_context_struct* ctx = (qmckl_context_struct*) context; 
+
+  qmckl_memory_struct* new_alloc = (qmckl_memory_struct*)
+    malloc(sizeof(qmckl_memory_struct));
+
+  if (new_alloc == NULL) {
+    qmckl_unlock(context);
+    return QMCKL_ALLOCATION_FAILED;
+  }
+
+  new_alloc->next    = NULL;
+  new_alloc->pointer = pointer;
+  new_alloc->size    = size;
+
+  qmckl_memory_struct* alloc = ctx->alloc;
+  if (alloc == NULL) {
+    ctx->alloc = new_alloc;
+  } else {
+    while (alloc != NULL) {
+      alloc = alloc->next;
+    }
+    alloc->next = new_alloc;
+  }
+
+  qmckl_unlock(context);
+
+  return QMCKL_SUCCESS;
+
+}
+
+
+
+
+ +
+

2.3 Remove memory

+
+

+The following function, called in qmckl_memory.c, removes a +pointer from the data structure. +It is forbidden to pass the NULL pointer. +If the context is QMCKL_NULL_CONTEXT, the function returns +immediately with QMCKL_SUCCESS. +

+ +
+
qmckl_exit_code qmckl_context_remove_memory(qmckl_context context,
+                                            const void* pointer) {
+  assert (pointer != NULL);
+
+  qmckl_lock(context);
+
+  if ( qmckl_context_check(context) == QMCKL_NULL_CONTEXT ) {
+    qmckl_unlock(context);
+    return QMCKL_SUCCESS;
+  }
+
+  qmckl_context_struct* ctx = (qmckl_context_struct*) context; 
+
+  qmckl_memory_struct* alloc = ctx->alloc;
+  qmckl_memory_struct* prev  = ctx->alloc;
+
+  while ( (alloc != NULL) && (alloc->pointer != pointer) ) {
+    prev = alloc;
+    alloc = alloc->next;
+  }
+
+  if (alloc != NULL) {
+    prev->next = alloc->next;
+    free(alloc);
+  } 
+
+  qmckl_unlock(context);
+
+  if (alloc != NULL) {
+    return QMCKL_SUCCESS;
+  } else {
+    return QMCKL_DEALLOCATION_FAILED;
+  }
+}
+
+
+
+
+
+ +
+

3 Error handling

+
+
+
+

3.1 Data structure

+
+
+
#define  QMCKL_MAX_FUN_LEN   256
+#define  QMCKL_MAX_MSG_LEN  1024
+
+typedef struct qmckl_error_struct {
+
+  qmckl_exit_code exit_code;
+  char function[QMCKL_MAX_FUN_LEN];
+  char message [QMCKL_MAX_MSG_LEN];
+
+} qmckl_error_struct;
+
+
+
+
+ +
+

3.2 Updating errors

+
+

+The error is updated in the context using +qmckl_context_update_error, although it is recommended to use +qmckl_context_set_error for the immutable variant. +When the error is set in the context, it is mandatory to specify +from which function the error is triggered, and a message +explaining the error. The exit code can't be QMCKL_SUCCESS. +

+ +
+
qmckl_exit_code
+qmckl_context_update_error(qmckl_context context,
+                           const qmckl_exit_code exit_code,
+                           const char* function_name,
+                           const char* message)
+{
+  /* Passing a function name and a message is mandatory. */
+  assert (function_name != NULL);
+  assert (message != NULL);
+
+  /* Exit codes are assumed valid. */
+  assert (exit_code >= 0);
+  assert (exit_code != QMCKL_SUCCESS);
+  assert (exit_code < QMCKL_INVALID_EXIT_CODE);
+
+  qmckl_lock(context);
+
+  /* The context is assumed to exist. */
+  assert (qmckl_context_check(context) != QMCKL_NULL_CONTEXT);
+
+  qmckl_context_struct* ctx = (qmckl_context_struct*) context;
+  assert (ctx != NULL); /* Impossible because the context is valid. */
+
+  if (ctx->error != NULL) {
+    free(ctx->error);
+    ctx->error = NULL;
+  }
+
+  qmckl_error_struct* error =
+    (qmckl_error_struct*) qmckl_malloc (context, sizeof(qmckl_error_struct));
+  error->exit_code = exit_code;
+  strcpy(error->function, function_name);
+  strcpy(error->message, message);
+
+  ctx->error = error;
+
+  qmckl_unlock(context);
+
+  return QMCKL_SUCCESS;
+}
+
+
+ +

+The qmckl_context_set_error function returns a new context with +the error domain updated. +

+ +
+
qmckl_context
+qmckl_context_set_error(qmckl_context context,
+                        const qmckl_exit_code exit_code,
+                        const char* function_name,
+                        const char* message)
+{
+  /* Passing a function name and a message is mandatory. */
+  assert (function_name != NULL);
+  assert (message != NULL);
+
+  /* Exit codes are assumed valid. */
+  assert (exit_code >= 0);
+  assert (exit_code != QMCKL_SUCCESS);
+  assert (exit_code < QMCKL_INVALID_EXIT_CODE);
+
+  /* The context is assumed to be valid */
+  if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT)
+    return QMCKL_NULL_CONTEXT;
+
+  qmckl_context new_context = qmckl_context_copy(context);
+
+  /* Should be impossible because the context is valid */
+  assert (new_context != QMCKL_NULL_CONTEXT);
+
+  if (qmckl_context_update_error(new_context,
+                                 exit_code,
+                                 function_name,
+                                 message) != QMCKL_SUCCESS) {
+    return context;
+  }
+
+  return new_context;
+}
+
+
+ + +

+To make a function fail, the qmckl_failwith function should be +called, such that information about the failure is stored in +the context. The desired exit code is given as an argument, as +well as the name of the function and an error message. The return +code of the function is the desired return code. +

+ +
+
qmckl_exit_code qmckl_failwith(qmckl_context context,
+                               const qmckl_exit_code exit_code,
+                               const char* function,
+                               const char* message) {
+
+  assert (exit_code > 0);
+  assert (exit_code < QMCKL_INVALID_EXIT_CODE);
+  assert (function != NULL);
+  assert (message  != NULL);
+  assert (strlen(function) < QMCKL_MAX_FUN_LEN);
+  assert (strlen(message)  < QMCKL_MAX_MSG_LEN);
+
+  if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT)
+    return QMCKL_NULL_CONTEXT;
+
+  const qmckl_exit_code rc = 
+    qmckl_context_update_error(context, exit_code, function, message);
+
+  assert (rc == QMCKL_SUCCESS);
+
+  return exit_code;
+}
+
+
+
+ +

+For example, this function can be used as +

+
+
if (x < 0) {
+  return qmckl_failwith(context,
+                        QMCKL_INVALID_ARG_2,
+                        "qmckl_function", 
+                        "Expected x >= 0");
+ }
+
+
+
+
+
+ +
+

4 Control of the numerical precision

+
+

+Controlling numerical precision enables optimizations. Here, the +default parameters determining the target numerical precision and +range are defined. +

+ + + + +++ ++ + + + + + + + + + + + +
QMCKL_DEFAULT_PRECISION53
QMCKL_DEFAULT_RANGE11
+ +
+
""" This script generates the C and Fortran constants for the error
+    codes from the org-mode table.
+"""
+
+result = [ "#+begin_src c :comments org :tangle (eval h)" ]
+for (text, code) in table:
+    text=text.replace("~","")
+    result += [ f"#define  {text:30s} {code:d}" ]
+result += [ "#+end_src" ]
+
+result += [ "" ]
+
+result += [ "#+begin_src f90 :comments org :tangle (eval fh) :exports none" ]
+for (text, code) in table:
+    text=text.replace("~","")
+    result += [ f"   integer, parameter :: {text:30s} = {code:d}" ]
+result += [ "#+end_src" ]
+
+return '\n'.join(result)
+
+
+
+ +
+
typedef struct qmckl_precision_struct {
+  int  precision;
+  int  range;
+} qmckl_precision_struct;
+
+
+ +

+The following functions set and get the required precision and +range. precision is an integer between 2 and 53, and range is an +integer between 2 and 11. +

+ +

+The setter functions functions return a new context as a 64-bit +integer. The getter functions return the value, as a 32-bit +integer. The update functions return QMCKL_SUCCESS or +QMCKL_FAILURE. +

+
+ +
+

4.1 Precision

+
+

+qmckl_context_update_precision modifies the parameter for the +numerical precision in a context. If the context doesn't have any +precision set yet, the default values are used. +

+ +
+
qmckl_exit_code qmckl_context_update_precision(const qmckl_context context, const int precision) {
+
+  if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT)
+    return QMCKL_INVALID_CONTEXT;
+
+if (precision <  2) {
+  return qmckl_failwith(context,
+                  QMCKL_INVALID_ARG_2,
+                  "qmckl_context_update_precision",
+                  "precision < 2");
+ }
+
+if (precision > 53) {
+  return qmckl_failwith(context,
+                  QMCKL_INVALID_ARG_2,
+                  "qmckl_context_update_precision",
+                  "precision > 53");
+ }
+
+qmckl_context_struct* ctx = (qmckl_context_struct*) context;
+
+/* This should be always true */
+assert (ctx != NULL);
+
+qmckl_lock(context);
+
+if (ctx->fp == NULL) {
+
+  ctx->fp = (qmckl_precision_struct*)
+    qmckl_malloc(context, sizeof(qmckl_precision_struct));
+
+  if (ctx->fp == NULL) {
+    return qmckl_failwith(context,
+                    QMCKL_ALLOCATION_FAILED,
+                    "qmckl_context_update_precision",
+                    "ctx->fp");
+  }
+
+  ctx->fp->precision = QMCKL_DEFAULT_PRECISION;
+  ctx->fp->range     = QMCKL_DEFAULT_RANGE;
+ }
+
+ctx->fp->precision = precision;
+
+qmckl_unlock(context);
+
+return QMCKL_SUCCESS;
+}
+
+
+ +
+
interface
+   integer (c_int32_t) function qmckl_context_update_precision(context, precision) bind(C)
+     use, intrinsic :: iso_c_binding
+     integer (c_int64_t), intent(in), value :: context
+     integer (c_int32_t), intent(in), value :: precision
+   end function qmckl_context_update_precision
+end interface
+
+
+ +

+qmckl_context_set_precision returns a copy of the context with a +different precision parameter. +

+ +
+
qmckl_context qmckl_context_set_precision(const qmckl_context context, const int precision) {
+  qmckl_context new_context = qmckl_context_copy(context);
+  if (new_context == 0) return 0;
+
+  if (qmckl_context_update_precision(new_context, precision) == QMCKL_FAILURE) return 0;
+
+  return new_context;
+}
+
+
+ +

+qmckl_context_get_precision returns the value of the numerical precision in the context. +

+ +
+
int qmckl_context_get_precision(const qmckl_context context) {
+  if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
+      return qmckl_failwith(context,
+                      QMCKL_INVALID_CONTEXT,
+                      "qmckl_context_get_precision",
+                      "");
+  }
+
+  const qmckl_context_struct* ctx = (qmckl_context_struct*) context;
+  if (ctx->fp != NULL) 
+    return ctx->fp->precision;
+  else
+    return QMCKL_DEFAULT_PRECISION;
+}
+
+
+ +
+
interface
+   integer (c_int32_t) function qmckl_context_get_precision(context) bind(C)
+     use, intrinsic :: iso_c_binding
+     integer (c_int64_t), intent(in), value :: context
+   end function qmckl_context_get_precision
+end interface
+
+
+
+
+ +
+

4.2 Range

+
+

+qmckl_context_update_range modifies the parameter for the numerical range in a given context. +

+ +
+
qmckl_exit_code qmckl_context_update_range(const qmckl_context context, const int range) {
+
+  if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT)
+    return QMCKL_INVALID_CONTEXT;
+
+  if (range <  2) {
+    return qmckl_failwith(context,
+                    QMCKL_INVALID_ARG_2,
+                    "qmckl_context_update_range",
+                    "range < 2");
+  }
+
+  if (range > 11) {
+    return qmckl_failwith(context,
+                    QMCKL_INVALID_ARG_2,
+                    "qmckl_context_update_range",
+                    "range > 11");
+  }
+
+  qmckl_context_struct* ctx = (qmckl_context_struct*) context;
+
+  /* This should be always true */
+  assert (ctx != NULL);
+
+  qmckl_lock(context);
+
+  if (ctx->fp == NULL) {
+
+    ctx->fp = (qmckl_precision_struct*)
+      qmckl_malloc(context, sizeof(qmckl_precision_struct));
+
+    if (ctx->fp == NULL) {
+      return qmckl_failwith(context,
+                      QMCKL_ALLOCATION_FAILED,
+                      "qmckl_context_update_range",
+                      "ctx->fp");
+    }
+
+    ctx->fp->precision = QMCKL_DEFAULT_PRECISION;
+    ctx->fp->range     = QMCKL_DEFAULT_RANGE;
+  }
+
+  ctx->fp->range = range;
+
+  qmckl_unlock(context);
+
+  return QMCKL_SUCCESS;
+}
+
+
+ +
+
interface
+   integer (c_int32_t) function qmckl_context_update_range(context, range) bind(C)
+     use, intrinsic :: iso_c_binding
+     integer (c_int64_t), intent(in), value :: context
+     integer (c_int32_t), intent(in), value :: range
+   end function qmckl_context_update_range
+end interface
+
+
+ +

+qmckl_context_set_range returns a copy of the context with a different precision parameter. +

+ +
+
qmckl_context qmckl_context_set_range(const qmckl_context context, const int range) {
+  qmckl_context new_context = qmckl_context_copy(context);
+  if (new_context == 0) return 0;
+
+  if (qmckl_context_update_range(new_context, range) == QMCKL_FAILURE) return 0;
+
+  return new_context;
+}
+
+
+ +

+qmckl_context_get_range returns the value of the numerical range in the context. +

+ +
+
int qmckl_context_get_range(const qmckl_context context) {
+  if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
+      return qmckl_failwith(context,
+                      QMCKL_INVALID_CONTEXT,
+                      "qmckl_context_get_range",
+                      "");
+  }
+
+  const qmckl_context_struct* ctx = (qmckl_context_struct*) context;
+  if (ctx->fp != NULL) 
+    return ctx->fp->range;
+  else
+    return QMCKL_DEFAULT_RANGE;
+}
+
+
+
+
+
+

4.3 Helper functions

+
+

+qmckl_context_get_epsilon returns \(\epsilon = 2^{1-n}\) where n is the precision. +

+ +
+
double qmckl_context_get_epsilon(const qmckl_context context) {
+  const int precision = qmckl_context_get_precision(context);
+  return 1. /  (double) (1L << (precision-1));
+}
+
+
+
+
+
+
+

5 TODO Basis set

+
+

+For H2 with the following basis set, +

+ +
+HYDROGEN
+S   5
+1         3.387000E+01           6.068000E-03
+2         5.095000E+00           4.530800E-02
+3         1.159000E+00           2.028220E-01
+4         3.258000E-01           5.039030E-01
+5         1.027000E-01           3.834210E-01
+S   1
+1         3.258000E-01           1.000000E+00
+S   1
+1         1.027000E-01           1.000000E+00
+P   1
+1         1.407000E+00           1.000000E+00
+P   1
+1         3.880000E-01           1.000000E+00
+D   1
+1         1.057000E+00           1.0000000
+
+ +

+we have: +

+ +
+type = 'G'
+shell_num = 12
+prim_num = 20
+SHELL_CENTER = [1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2]
+SHELL_ANG_MOM = ['S', 'S', 'S', 'P', 'P', 'D', 'S', 'S', 'S', 'P', 'P', 'D']
+SHELL_PRIM_NUM = [5, 1, 1, 1, 1, 1, 5, 1, 1, 1, 1, 1]
+prim_index = [1, 6, 7, 8, 9, 10, 11, 16, 17, 18, 19, 20]
+EXPONENT = [ 33.87, 5.095, 1.159, 0.3258, 0.1027, 0.3258, 0.1027,
+             1.407, 0.388, 1.057, 33.87, 5.095, 1.159, 0.3258, 0.1027,
+             0.3258, 0.1027, 1.407, 0.388, 1.057]
+COEFFICIENT = [ 0.006068, 0.045308, 0.202822, 0.503903, 0.383421,
+                1.0, 1.0, 1.0, 1.0, 1.0, 0.006068, 0.045308, 0.202822,
+                0.503903, 0.383421, 1.0, 1.0, 1.0, 1.0, 1.0]
+
+
+ +
+

5.1 Data structure

+
+
+
typedef struct qmckl_ao_basis_struct {
+
+  int64_t   shell_num;
+  int64_t   prim_num;
+  int64_t * shell_center;
+  int32_t * shell_ang_mom;
+  double  * shell_factor;
+  double  * exponent    ;
+  double  * coefficient ;
+  int64_t * shell_prim_num;
+  char      type;
+
+} qmckl_ao_basis_struct;
+
+
+
+
+ +
+

5.2 qmckl_context_update_ao_basis

+
+

+Updates the data describing the AO basis set into the context. +

+ + + + +++ ++ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
typeGaussian or Slater
shell_numNumber of shells
prim_numTotal number of primitives
SHELL_CENTER(shell_num)Id of the nucleus on which the shell is centered
SHELL_ANG_MOM(shell_num)Id of the nucleus on which the shell is centered
SHELL_FACTOR(shell_num)Normalization factor for the shell
SHELL_PRIM_NUM(shell_num)Number of primitives in the shell
SHELL_PRIM_INDEX(shell_num)Address of the first primitive of the shelll in the EXPONENT array
EXPONENT(prim_num)Array of exponents
COEFFICIENT(prim_num)Array of coefficients
+ +
+
qmckl_exit_code
+qmckl_context_update_ao_basis(qmckl_context   context         ,
+                              const char      type            ,
+                              const int64_t   shell_num       ,
+                              const int64_t   prim_num        ,
+                              const int64_t * SHELL_CENTER    ,
+                              const int32_t * SHELL_ANG_MOM   ,
+                              const double  * SHELL_FACTOR    ,
+                              const int64_t * SHELL_PRIM_NUM  ,
+                              const int64_t * SHELL_PRIM_INDEX,
+                              const double  * EXPONENT        ,
+                              const double  * COEFFICIENT);
+
+
+
+ +
+

5.2.1 Source

+
+
+
qmckl_exit_code
+qmckl_context_update_ao_basis(qmckl_context   context     , const char      type,
+                              const int64_t   shell_num   , const int64_t   prim_num,
+                              const int64_t * SHELL_CENTER, const int32_t * SHELL_ANG_MOM,
+                              const double  * SHELL_FACTOR, const int64_t * SHELL_PRIM_NUM,
+                              const int64_t * SHELL_PRIM_INDEX,
+                              const double  * EXPONENT    , const double  * COEFFICIENT)
+{
+
+  int64_t i;
+
+  /* Check input */
+
+  if (type != 'G' && type != 'S') return QMCKL_FAILURE;
+  if (shell_num <= 0) return QMCKL_FAILURE;
+  if (prim_num  <= 0) return QMCKL_FAILURE;
+  if (prim_num  < shell_num) return QMCKL_FAILURE;
+
+  for (i=0 ; i<shell_num ; i++) {
+    if (SHELL_CENTER[i] <= 0) return QMCKL_FAILURE;
+    if (SHELL_PRIM_NUM[i] <= 0) return QMCKL_FAILURE;
+    if (SHELL_ANG_MOM[i] < 0) return QMCKL_FAILURE;
+    if (SHELL_PRIM_INDEX[i] < 0) return QMCKL_FAILURE;
+  }
+
+  for (i=0 ; i<prim_num ; i++) {
+    if (EXPONENT[i] <= 0) return QMCKL_FAILURE;
+  }
+
+  qmckl_context_struct* ctx = (qmckl_context_struct*) context;
+  if (ctx == NULL) return QMCKL_FAILURE;
+
+  qmckl_ao_basis_struct* basis =
+    (qmckl_ao_basis_struct*) qmckl_malloc (context, sizeof(qmckl_ao_basis_struct));
+  if (basis == NULL) return QMCKL_ALLOCATION_FAILED;
+
+
+  /* Memory allocations */
+
+  basis->shell_center  = (int64_t*) qmckl_malloc (context, shell_num * sizeof(int64_t));
+  if (basis->shell_center == NULL) {
+    qmckl_free(context, basis);
+    basis = NULL;
+    return QMCKL_FAILURE;
+  }
+
+  basis->shell_ang_mom = (int32_t*) qmckl_malloc (context, shell_num * sizeof(int32_t));
+  if (basis->shell_ang_mom == NULL) {
+    qmckl_free(context, basis->shell_center);
+    basis->shell_center = NULL;
+    qmckl_free(context, basis);
+    basis = NULL;
+    return QMCKL_FAILURE;
+  }
+
+  basis->shell_prim_num= (int64_t*) qmckl_malloc (context, shell_num * sizeof(int64_t));
+  if (basis->shell_prim_num == NULL) {
+    qmckl_free(context, basis->shell_ang_mom);
+    basis->shell_ang_mom = NULL;
+    qmckl_free(context, basis->shell_center);
+    basis->shell_center = NULL;
+    qmckl_free(context, basis);
+    basis = NULL;
+    return QMCKL_FAILURE;
+  }
+
+  basis->shell_factor  = (double *) qmckl_malloc (context, shell_num * sizeof(double ));
+  if (basis->shell_factor == NULL) {
+    qmckl_free(context, basis->shell_prim_num);
+    basis->shell_prim_num = NULL;
+    qmckl_free(context, basis->shell_ang_mom);
+    basis->shell_ang_mom = NULL;
+    qmckl_free(context, basis->shell_center);
+    basis->shell_center = NULL;
+    qmckl_free(context, basis);
+    basis = NULL;
+    return QMCKL_FAILURE;
+  }
+
+  basis->exponent      = (double *) qmckl_malloc (context, prim_num  * sizeof(double ));
+  if (basis->exponent == NULL) {
+    qmckl_free(context, basis->shell_factor);
+    basis->shell_factor = NULL;
+    qmckl_free(context, basis->shell_prim_num);
+    basis->shell_prim_num = NULL;
+    qmckl_free(context, basis->shell_ang_mom);
+    basis->shell_ang_mom = NULL;
+    qmckl_free(context, basis->shell_center);
+    basis->shell_center = NULL;
+    qmckl_free(context, basis);
+    basis = NULL;
+    return QMCKL_FAILURE;
+  }
+
+  basis->coefficient   = (double *) qmckl_malloc (context, prim_num  * sizeof(double ));
+  if (basis->coefficient == NULL) {
+    qmckl_free(context, basis->exponent);
+    basis->exponent = NULL;
+    qmckl_free(context, basis->shell_factor);
+    basis->shell_factor = NULL;
+    qmckl_free(context, basis->shell_prim_num);
+    basis->shell_prim_num = NULL;
+    qmckl_free(context, basis->shell_ang_mom);
+    basis->shell_ang_mom = NULL;
+    qmckl_free(context, basis->shell_center);
+    basis->shell_center = NULL;
+    qmckl_free(context, basis);
+    basis = NULL;
+    return QMCKL_FAILURE;
+  }
+
+
+  /* Assign data */
+
+  basis->type      = type;
+  basis->shell_num = shell_num;
+  basis->prim_num  = prim_num;
+
+  for (i=0 ; i<shell_num ; i++) {
+    basis->shell_center  [i] = SHELL_CENTER  [i];
+    basis->shell_ang_mom [i] = SHELL_ANG_MOM [i];
+    basis->shell_prim_num[i] = SHELL_PRIM_NUM[i];
+    basis->shell_factor  [i] = SHELL_FACTOR  [i];
+  }
+
+  for (i=0 ; i<prim_num ; i++) {
+    basis->exponent   [i] = EXPONENT[i];
+    basis->coefficient[i] = COEFFICIENT[i];
+  }
+
+  ctx->ao_basis = basis;
+  return QMCKL_SUCCESS;
+}
+
+
+
+
+ +
+

5.2.2 Fortran interface

+
+
+
interface
+   integer (c_int32_t) function qmckl_context_update_ao_basis(context, &
+        typ, shell_num, prim_num, SHELL_CENTER, SHELL_ANG_MOM, SHELL_FACTOR, &
+        SHELL_PRIM_NUM, SHELL_PRIM_INDEX, EXPONENT, COEFFICIENT) bind(C)
+     use, intrinsic :: iso_c_binding
+     integer (c_int64_t), intent(in), value :: context
+     character(c_char)  , intent(in), value :: typ
+     integer (c_int64_t), intent(in), value :: shell_num
+     integer (c_int64_t), intent(in), value :: prim_num
+     integer (c_int64_t), intent(in)        :: SHELL_CENTER(shell_num)
+     integer (c_int32_t), intent(in)        :: SHELL_ANG_MOM(shell_num)
+     double precision   , intent(in)        :: SHELL_FACTOR(shell_num)
+     integer (c_int64_t), intent(in)        :: SHELL_PRIM_NUM(shell_num)
+     integer (c_int64_t), intent(in)        :: SHELL_PRIM_INDEX(shell_num)
+     double precision   , intent(in)        :: EXPONENT(prim_num)
+     double precision   , intent(in)        :: COEFFICIENT(prim_num)
+   end function qmckl_context_update_ao_basis
+end interface
+
+
+
+
+ +
+

5.2.3 TODO Test

+
+
+ +
+

5.3 qmckl_context_set_ao_basis

+
+

+Sets the data describing the AO basis set into the context. +

+ + + + +++ ++ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
typeGaussian or Slater
shell_numNumber of shells
prim_numTotal number of primitives
SHELL_CENTER(shell_num)Id of the nucleus on which the shell is centered
SHELL_ANG_MOM(shell_num)Id of the nucleus on which the shell is centered
SHELL_FACTOR(shell_num)Normalization factor for the shell
SHELL_PRIM_NUM(shell_num)Number of primitives in the shell
SHELL_PRIM_INDEX(shell_num)Address of the first primitive of the shelll in the EXPONENT array
EXPONENT(prim_num)Array of exponents
COEFFICIENT(prim_num)Array of coefficients
+ +
+
qmckl_context
+qmckl_context_set_ao_basis(const qmckl_context context     , const char      type,
+                           const int64_t       shell_num   , const int64_t   prim_num,
+                           const int64_t *     SHELL_CENTER, const int32_t * SHELL_ANG_MOM,
+                           const double  *     SHELL_FACTOR, const int64_t * SHELL_PRIM_NUM,
+                           const int64_t *     SHELL_PRIM_INDEX,
+                           const double  *     EXPONENT    , const double  * COEFFICIENT);
+
+
+
+ +
+

5.3.1 Source

+
+
+
qmckl_context
+qmckl_context_set_ao_basis(const qmckl_context context     , const char      type,
+                           const int64_t       shell_num   , const int64_t   prim_num,
+                           const int64_t *     SHELL_CENTER, const int32_t * SHELL_ANG_MOM,
+                           const double  *     SHELL_FACTOR, const int64_t * SHELL_PRIM_NUM,
+                           const int64_t *     SHELL_PRIM_INDEX,
+                           const double  *     EXPONENT    , const double  * COEFFICIENT)
+{
+
+  qmckl_context new_context = qmckl_context_copy(context);
+  if (new_context == 0) return 0;
+
+  if (qmckl_context_update_ao_basis(new_context, type, shell_num, prim_num,
+                                    SHELL_CENTER, SHELL_ANG_MOM, SHELL_FACTOR,
+                                    SHELL_PRIM_NUM, SHELL_PRIM_INDEX, EXPONENT,
+                                    COEFFICIENT
+                                    ) == QMCKL_FAILURE)
+    return 0;
+
+  return new_context;
+}
+
+
+
+
+ +
+

5.3.2 Fortran interface

+
+
+
interface
+   integer (c_int64_t) function qmckl_context_set_ao_basis(context, &
+        typ, shell_num, prim_num, SHELL_CENTER, SHELL_ANG_MOM, SHELL_FACTOR, &
+        SHELL_PRIM_NUM, SHELL_PRIM_INDEX, EXPONENT, COEFFICIENT) bind(C)
+     use, intrinsic :: iso_c_binding
+     integer (c_int64_t), intent(in), value :: context
+     character(c_char)  , intent(in), value :: typ
+     integer (c_int64_t), intent(in), value :: shell_num
+     integer (c_int64_t), intent(in), value :: prim_num
+     integer (c_int64_t), intent(in)        :: SHELL_CENTER(shell_num)
+     integer (c_int32_t), intent(in)        :: SHELL_ANG_MOM(shell_num)
+     double precision   , intent(in)        :: SHELL_FACTOR(shell_num)
+     integer (c_int64_t), intent(in)        :: SHELL_PRIM_NUM(shell_num)
+     integer (c_int64_t), intent(in)        :: SHELL_PRIM_INDEX(shell_num)
+     double precision   , intent(in)        :: EXPONENT(prim_num)
+     double precision   , intent(in)        :: COEFFICIENT(prim_num)
+   end function qmckl_context_set_ao_basis
+end interface
+
+
+
+
+ +
+

5.3.3 TODO Test

+
+
+
+
+
+

Author: TREX CoE

+

Created: 2021-03-19 Fri 12:49

+

Validate

+
+ + diff --git a/qmckl_distance.html b/qmckl_distance.html new file mode 100644 index 0000000..0798742 --- /dev/null +++ b/qmckl_distance.html @@ -0,0 +1,628 @@ + + + + + + + +Distances + + + + + + + + + + + + + +
+ UP + | + HOME +
+

Distances

+
+

Table of Contents

+
+ +
+
+

+Functions for the computation of distances between particles. +

+ +
+

1 Squared distance

+
+

+qmckl_distance_sq computes the matrix of the squared distances + between all pairs of points in two sets, one point within each set: +

+ +

+\[ + C_{ij} = \sum_{k=1}^3 (A_{k,i}-B_{k,j})^2 + \] +

+ + + + +++ ++ ++ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
contextinputGlobal state
transainputArray A is N: Normal, T: Transposed
transbinputArray B is N: Normal, T: Transposed
minputNumber of points in the first set
ninputNumber of points in the second set
A(lda,3)inputArray containing the \(m \times 3\) matrix \(A\)
ldainputLeading dimension of array A
B(ldb,3)inputArray containing the \(n \times 3\) matrix \(B\)
ldbinputLeading dimension of array B
C(ldc,n)outputArray containing the \(m \times n\) matrix \(C\)
ldcinputLeading dimension of array C
+
+ +
+

1.0.1 Requirements

+
+
    +
  • context is not 0
    • +
  • m > 0
    • +
  • n > 0
    • +
  • lda >= 3 if transa is N
    • +
  • lda >= m if transa is T
    • +
  • ldb >= 3 if transb is N
    • +
  • ldb >= n if transb is T
    • +
  • ldc >= m
    • +
  • A is allocated with at least \(3 \times m \times 8\) bytes
    • +
  • B is allocated with at least \(3 \times n \times 8\) bytes
    • +
  • C is allocated with at least \(m \times n \times 8\) bytes
    • +
+
+
+ +
+

1.0.2 Performance

+
+

+This function might be more efficient when A and B are +transposed. +

+ +
+
qmckl_exit_code qmckl_distance_sq(const qmckl_context context,
+                                  const char transa, const char transb,
+                                  const int64_t m, const int64_t n,
+                                  const double *A, const int64_t lda,
+                                  const double *B, const int64_t ldb,
+                                  const double *C, const int64_t ldc);
+
+
+
+
+ +
+

1.0.3 Source

+
+
+
integer function qmckl_distance_sq_f(context, transa, transb, m, n, A, LDA, B, LDB, C, LDC) result(info)
+  use qmckl
+  implicit none
+  integer*8  , intent(in)  :: context
+  character  , intent(in)  :: transa, transb
+  integer*8  , intent(in)  :: m, n
+  integer*8  , intent(in)  :: lda
+  real*8     , intent(in)  :: A(lda,*)
+  integer*8  , intent(in)  :: ldb
+  real*8     , intent(in)  :: B(ldb,*)
+  integer*8  , intent(in)  :: ldc
+  real*8     , intent(out) :: C(ldc,*)
+
+  integer*8 :: i,j
+  real*8    :: x, y, z
+  integer   :: transab
+
+  info = 0
+
+  if (context == QMCKL_NULL_CONTEXT) then
+     info = QMCKL_INVALID_CONTEXT
+     return
+  endif
+
+  if (m <= 0_8) then
+     info = QMCKL_INVALID_ARG_4
+     return
+  endif
+
+  if (n <= 0_8) then
+     info = QMCKL_INVALID_ARG_5
+     return
+  endif
+
+  if (transa == 'N' .or. transa == 'n') then
+     transab = 0
+  else if (transa == 'T' .or. transa == 't') then
+     transab = 1
+  else
+     transab = -100
+  endif
+
+  if (transb == 'N' .or. transb == 'n') then
+     continue
+  else if (transa == 'T' .or. transa == 't') then
+     transab = transab + 2
+  else
+     transab = -100
+  endif
+
+  if (transab < 0) then
+     info = QMCKL_INVALID_ARG_1
+     return 
+  endif
+
+  if (iand(transab,1) == 0 .and. LDA < 3) then
+     info = QMCKL_INVALID_ARG_7
+     return
+  endif
+
+  if (iand(transab,1) == 1 .and. LDA < m) then
+     info = QMCKL_INVALID_ARG_7
+     return
+  endif
+
+  if (iand(transab,2) == 0 .and. LDA < 3) then
+     info = QMCKL_INVALID_ARG_7
+     return
+  endif
+
+  if (iand(transab,2) == 2 .and. LDA < m) then
+     info = QMCKL_INVALID_ARG_7
+     return
+  endif
+
+
+  select case (transab)
+
+  case(0)
+
+     do j=1,n
+        do i=1,m
+           x = A(1,i) - B(1,j)
+           y = A(2,i) - B(2,j)
+           z = A(3,i) - B(3,j)
+           C(i,j) = x*x + y*y + z*z
+        end do
+     end do
+
+  case(1)
+
+     do j=1,n
+        do i=1,m
+           x = A(i,1) - B(1,j)
+           y = A(i,2) - B(2,j)
+           z = A(i,3) - B(3,j)
+           C(i,j) = x*x + y*y + z*z
+        end do
+     end do
+
+  case(2)
+
+     do j=1,n
+        do i=1,m
+           x = A(1,i) - B(j,1)
+           y = A(2,i) - B(j,2)
+           z = A(3,i) - B(j,3)
+           C(i,j) = x*x + y*y + z*z
+        end do
+     end do
+
+  case(3)
+
+     do j=1,n
+        do i=1,m
+           x = A(i,1) - B(j,1)
+           y = A(i,2) - B(j,2)
+           z = A(i,3) - B(j,3)
+           C(i,j) = x*x + y*y + z*z
+        end do
+     end do
+
+  end select
+
+end function qmckl_distance_sq_f
+
+
+
+
+
+
+
+

Author: TREX CoE

+

Created: 2021-03-19 Fri 12:49

+

Validate

+
+ + diff --git a/qmckl_error.html b/qmckl_error.html new file mode 100644 index 0000000..bbc71ea --- /dev/null +++ b/qmckl_error.html @@ -0,0 +1,452 @@ + + + + + + + +Error handling + + + + + + + + + + + +
+ UP + | + HOME +
+

Error handling

+
+

Table of Contents

+
+
    +
    • +
+
+
+ +
+

+
+

+The library should never make the calling programs abort, nor +perform any input/output operations. This decision has to be taken +by the developer of the code calling the library. +

+ +

+All the functions return with an exit code, defined as +

+
+
typedef int32_t qmckl_exit_code;
+
+
+ +

+The exit code returns the completion status of the function to the +calling program. When a function call completed successfully, +QMCKL_SUCCESS is returned. If one of the functions of +the library fails to complete the requested task, an appropriate +error code is returned to the program. +

+ +

+Here is the complete list of exit codes. +

+ + + + +++ ++ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
QMCKL_SUCCESS0
QMCKL_INVALID_ARG_11
QMCKL_INVALID_ARG_22
QMCKL_INVALID_ARG_33
QMCKL_INVALID_ARG_44
QMCKL_INVALID_ARG_55
QMCKL_INVALID_ARG_66
QMCKL_INVALID_ARG_77
QMCKL_INVALID_ARG_88
QMCKL_INVALID_ARG_99
QMCKL_INVALID_ARG_1010
QMCKL_FAILURE101
QMCKL_ERRNO102
QMCKL_INVALID_CONTEXT103
QMCKL_ALLOCATION_FAILED104
QMCKL_DEALLOCATION_FAILED105
QMCKL_INVALID_EXIT_CODE106
+
+
+
+
+

Author: TREX CoE

+

Created: 2021-03-19 Fri 12:49

+

Validate

+
+ + diff --git a/qmckl_memory.html b/qmckl_memory.html new file mode 100644 index 0000000..ab80e0d --- /dev/null +++ b/qmckl_memory.html @@ -0,0 +1,436 @@ + + + + + + + +Memory management + + + + + + + + + + + +
+ UP + | + HOME +
+

Memory management

+
+

Table of Contents

+
+
    +
    • +
+
+
+

+We override the allocation functions to enable the possibility of +optimized libraries to fine-tune the memory allocation. +

+ + +
+

+
+

+Memory allocation inside the library should be done with +qmckl_malloc. It lets the library choose how the memory will be +allocated, and a pointer is returned to the user. The context is +passed to let the library store data related to the allocation +inside the context. In this particular implementation of the library, +we store a list of allocated pointers so that all the memory can be +properly freed when the library is de-initialized. +If the allocation failed, the NULL pointer is returned. +

+ +
+
void* qmckl_malloc(qmckl_context context,
+                   const size_t size);
+
+
+ +
+
void* qmckl_malloc(qmckl_context context, const size_t size) {
+
+  void * pointer = malloc( (size_t) size );
+
+  if (qmckl_context_check(context) != QMCKL_NULL_CONTEXT) {
+    qmckl_exit_code rc;
+    rc = qmckl_context_append_memory(context, pointer, size);
+    assert (rc == QMCKL_SUCCESS);
+  }
+
+  return pointer;
+}
+
+
+ +
+
interface
+   type (c_ptr) function qmckl_malloc (context, size) bind(C)
+     use, intrinsic :: iso_c_binding
+     integer (c_int64_t), intent(in), value :: context
+     integer (c_int64_t), intent(in), value :: size
+   end function qmckl_malloc
+end interface
+
+
+ +
+
qmckl_context context = qmckl_context_create();    
+
+int *a = (int*) qmckl_malloc(context, 3*sizeof(int));
+munit_assert(a != NULL);
+
+a[0] = 1;  munit_assert_int(a[0], ==, 1);
+a[1] = 2;  munit_assert_int(a[1], ==, 2);
+a[2] = 3;  munit_assert_int(a[2], ==, 3);
+
+
+ +

+When freeing the memory with qmckl_free, the context is passed, in +case some important information has been stored related to memory +allocation and needs to be updated. +

+ +
+
qmckl_exit_code qmckl_free(qmckl_context context,
+                           void *ptr);
+
+
+ +
+
interface
+   integer (c_int32_t) function qmckl_free (context, ptr) bind(C)
+     use, intrinsic :: iso_c_binding
+     integer (c_int64_t), intent(in), value :: context
+     type (c_ptr), intent(in), value :: ptr
+   end function qmckl_free
+end interface
+
+
+ +
+
qmckl_exit_code qmckl_free(qmckl_context context, void *ptr) {
+  if (qmckl_context_check(context) != QMCKL_NULL_CONTEXT) {
+
+    if (ptr == NULL) {
+      return qmckl_failwith(context,
+                            QMCKL_INVALID_ARG_2,
+                            "qmckl_free",
+                            "NULL pointer");
+    }
+
+    qmckl_exit_code rc;
+    rc = qmckl_context_remove_memory(context, ptr);
+
+    assert (rc == QMCKL_SUCCESS);
+  }
+  free(ptr);
+  return QMCKL_SUCCESS;
+}
+
+
+
+
+
+
+

Author: TREX CoE

+

Created: 2021-03-19 Fri 12:49

+

Validate

+
+ + diff --git a/test_qmckl.html b/test_qmckl.html new file mode 100644 index 0000000..1709f85 --- /dev/null +++ b/test_qmckl.html @@ -0,0 +1,318 @@ + + + + + + + +Testing + + + + + + + + + + + +
+ UP + | + HOME +
+

Testing

+ +
+
+

Author: TREX CoE

+

Created: 2021-03-19 Fri 12:49

+

Validate

+
+ + diff --git a/theme.setup b/theme.setup new file mode 100644 index 0000000..da0d396 --- /dev/null +++ b/theme.setup @@ -0,0 +1,15 @@ +# -*- mode: org; -*- + +#+HTML_LINK_HOME: index.html +#+OPTIONS: H:4 num:t toc:t \n:nil @:t ::t |:t ^:t -:t f:t *:t <:t d:(HIDE) + +# SETUPFILE: ../docs/org-html-themes/org/theme-readtheorg.setup + +#+INFOJS_OPT: toc:t mouse:underline path:org-info.js +#+HTML_HEAD: + +#+STARTUP: align fold nodlcheck hidestars oddeven lognotestate +#+AUTHOR: TREX CoE +#+LANGUAGE: en + +