Fixed address sanitizer

org/qmckl_context.org

@@ -222,6 +222,7 @@ qmckl_context_touch(const qmckl_context context)
 
   qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
 
+  ctx->electron.walker_old = ctx->electron.walker;
   ctx->date += 1UL;
   ctx->point.date = ctx-> date;
   return QMCKL_SUCCESS;

org/qmckl_jastrow_champ.org

@@ -471,7 +471,6 @@ if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
 qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
 
 if (mask != 0 && !(ctx->jastrow_champ.uninitialized & mask)) {
-  printf("%d %d\n", mask, ctx->jastrow_champ.uninitialized );
     return qmckl_failwith( context,
                            QMCKL_ALREADY_SET,
                            "qmckl_set_jastrow_champ_*",
@@ -928,10 +927,6 @@ qmckl_exit_code qmckl_finalize_jastrow_champ(qmckl_context context) {
  /* Check for the necessary information */
  /* ----------------------------------- */
 
-  /* Check for the electron data
-     1. elec_num
-     2. ee_distances_rescaled
-  ,*/
   if (!(ctx->electron.provided)) {
     return qmckl_failwith( context,
                            QMCKL_NOT_PROVIDED,
@@ -939,10 +934,6 @@ qmckl_exit_code qmckl_finalize_jastrow_champ(qmckl_context context) {
                            NULL);
   }
 
-  /* Check for the nucleus data
-     1. nucl_num
-     2. en_distances_rescaled
-  ,*/
   if (!(ctx->nucleus.provided)) {
     return qmckl_failwith( context,
                            QMCKL_NOT_PROVIDED,
@@ -4172,6 +4163,11 @@ integer function qmckl_compute_jastrow_champ_factor_en_doc_f( &
      return
   endif
 
+  if (type_nucl_num <= 0) then
+     info = QMCKL_INVALID_ARG_4
+     return
+  endif
+
   if (aord_num < 0) then
      info = QMCKL_INVALID_ARG_7
      return
@@ -4249,94 +4245,6 @@ end function qmckl_compute_jastrow_champ_factor_en_doc_f
    end function qmckl_compute_jastrow_champ_factor_en_doc
     #+end_src
 
-     #+begin_src c :comments org :tangle (eval c) :noweb yes
-/*
-qmckl_exit_code qmckl_compute_jastrow_champ_factor_en (
-      const qmckl_context context,
-      const int64_t walk_num,
-      const int64_t elec_num,
-      const int64_t nucl_num,
-      const int64_t type_nucl_num,
-      const int64_t* type_nucl_vector,
-      const int64_t aord_num,
-      const double* a_vector,
-      const double* en_distance_rescaled,
-      const double* asymp_jasa,
-      double* const factor_en ) {
-
-  double  x, x1, power_ser;
-
-
-  if (context == QMCKL_NULL_CONTEXT) {
-     return QMCKL_INVALID_CONTEXT;
-  }
-
-  if (walk_num <= 0) {
-     return QMCKL_INVALID_ARG_2;
-  }
-
-  if (elec_num <= 0) {
-     return QMCKL_INVALID_ARG_3;
-  }
-
-  if (nucl_num <= 0) {
-     return QMCKL_INVALID_ARG_4;
-  }
-
-  if (type_nucl_num <= 0) {
-    return QMCKL_INVALID_ARG_5;
-  }
-
-  if (type_nucl_vector == NULL) {
-    return QMCKL_INVALID_ARG_6;
-  }
-
-  if (aord_num < 0) {
-     return QMCKL_INVALID_ARG_7;
-  }
-
-  if (a_vector == NULL) {
-    return QMCKL_INVALID_ARG_8;
-  }
-
-  if (en_distance_rescaled == NULL) {
-    return QMCKL_INVALID_ARG_9;
-  }
-
-  if (factor_en == NULL) {
-    return QMCKL_INVALID_ARG_10;
-  }
-
-
-  for (int nw = 0; nw < walk_num; ++nw ) {
-    // init array
-    factor_en[nw] = 0.0;
-    for (int a = 0; a < nucl_num; ++a ) {
-      for (int i = 0; i < elec_num; ++i ) {
-        x = en_distance_rescaled[i + a * elec_num + nw * (elec_num * nucl_num)];
-        x1 = x;
-        power_ser = 0.0;
-
-        for (int p = 2; p < aord_num+1; ++p) {
-          x = x * x1;
-          power_ser = power_ser + a_vector[p+ (type_nucl_vector[a]-1) * aord_num] * x;
-        }
-
-        factor_en[nw] = factor_en[nw] + a_vector[0 + (type_nucl_vector[a]-1)*aord_num] * x1 / \
-                        (1.0 + a_vector[1 + (type_nucl_vector[a]-1) * aord_num] * x1) + \
-                        power_ser;
-
-      }
-      factor_en[nw] = factor_en[nw] + asymp_jasa[type_nucl_vector[a];
-    }
-  }
-
-  return QMCKL_SUCCESS;
-}
-*/
-     #+end_src
-
-
     #+CALL: generate_c_header(table=qmckl_factor_en_args,rettyp=get_value("CRetType"),fname=get_value("Name"))
 
     #+RESULTS:
@@ -5802,7 +5710,9 @@ qmckl_exit_code qmckl_compute_een_rescaled_e_hpc (
   // probably in C is better [cord+1, Ne*(Ne-1)/2]
   // elec_pairs = (elec_num * (elec_num - 1)) / 2;
   // len_een_ij = elec_pairs * (cord_num + 1);
-  double* een_rescaled_e_ij = (double *) malloc (len_een_ij * sizeof(double));
+  qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
+  mem_info.size = len_een_ij * sizeof(double);
+  double* een_rescaled_e_ij = (double*) qmckl_malloc(context, mem_info);
   if (een_rescaled_e_ij == NULL) {
     return qmckl_failwith( context,
                            QMCKL_ALLOCATION_FAILED,
@@ -5820,7 +5730,7 @@ qmckl_exit_code qmckl_compute_een_rescaled_e_hpc (
       }
     }
   }
-  */
+  ,*/
 
   for (int nw = 0; nw < walk_num; ++nw) {
 
@@ -5878,9 +5788,9 @@ qmckl_exit_code qmckl_compute_een_rescaled_e_hpc (
 
   }
 
-  free(een_rescaled_e_ij);
+  qmckl_exit_code rc = qmckl_free(context,een_rescaled_e_ij);
 
-  return QMCKL_SUCCESS;
+  return rc;
 }
      #+end_src
 
@@ -7262,7 +7172,7 @@ qmckl_exit_code qmckl_get_jastrow_champ_tmp_c(qmckl_context context, double* con
   assert (ctx != NULL);
 
   size_t sze = (ctx->jastrow_champ.cord_num) * (ctx->jastrow_champ.cord_num + 1)
-               * ctx->electron.num * ctx->nucleus.num * ctx->electron.walker.num;
+               ,* ctx->electron.num * ctx->nucleus.num * ctx->electron.walker.num;
   memcpy(tmp_c, ctx->jastrow_champ.tmp_c, sze * sizeof(double));
 
   return QMCKL_SUCCESS;
@@ -7286,7 +7196,7 @@ qmckl_exit_code qmckl_get_jastrow_champ_dtmp_c(qmckl_context context, double* co
   assert (ctx != NULL);
 
   size_t sze = (ctx->jastrow_champ.cord_num) * (ctx->jastrow_champ.cord_num + 1)
-               *4* ctx->electron.num * ctx->nucleus.num * ctx->electron.walker.num;
+    ,*4* ctx->electron.num * ctx->nucleus.num * ctx->electron.walker.num;
   memcpy(dtmp_c, ctx->jastrow_champ.dtmp_c, sze * sizeof(double));
 
   return QMCKL_SUCCESS;
@@ -7313,9 +7223,23 @@ qmckl_exit_code qmckl_provide_jastrow_champ_c_vector_full(qmckl_context context)
   qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
   assert (ctx != NULL);
 
+  qmckl_exit_code rc = QMCKL_SUCCESS;
+  
  /* Compute if necessary */
   if (ctx->date > ctx->jastrow_champ.c_vector_full_date) {
 
+    if (ctx->electron.walker.num > ctx->electron.walker_old.num) {
+      if (ctx->jastrow_champ.c_vector_full != NULL) {
+        rc = qmckl_free(context, ctx->jastrow_champ.c_vector_full);
+        if (rc != QMCKL_SUCCESS) {
+          return qmckl_failwith( context, rc,
+                                 "qmckl_provide_jastrow_champ_c_vector_full",
+                                 "Unable to free ctx->jastrow_champ.c_vector_full");
+        }
+        ctx->jastrow_champ.c_vector_full = NULL;
+      }
+    }
+
  /* Allocate array */
     if (ctx->jastrow_champ.c_vector_full == NULL) {
 
@@ -7332,7 +7256,6 @@ qmckl_exit_code qmckl_provide_jastrow_champ_c_vector_full(qmckl_context context)
       ctx->jastrow_champ.c_vector_full = c_vector_full;
     }
 
-    qmckl_exit_code rc;
     rc = qmckl_compute_c_vector_full(context,
                                       ctx->nucleus.num,
                                       ctx->jastrow_champ.dim_c_vector,
@@ -7360,9 +7283,23 @@ qmckl_exit_code qmckl_provide_lkpm_combined_index(qmckl_context context)
   qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
   assert (ctx != NULL);
 
+  qmckl_exit_code rc = QMCKL_SUCCESS;
+
  /* Compute if necessary */
   if (ctx->date > ctx->jastrow_champ.lkpm_combined_index_date) {
 
+    if (ctx->electron.walker.num > ctx->electron.walker_old.num) {
+      if (ctx->jastrow_champ.lkpm_combined_index != NULL) {
+        rc = qmckl_free(context, ctx->jastrow_champ.lkpm_combined_index);
+        if (rc != QMCKL_SUCCESS) {
+          return qmckl_failwith( context, rc,
+                                 "qmckl_provide_jastrow_champ_factor_ee",
+                                 "Unable to free ctx->jastrow_champ.lkpm_combined_index");
+        }
+        ctx->jastrow_champ.lkpm_combined_index = NULL;
+      }
+    }
+
  /* Allocate array */
     if (ctx->jastrow_champ.lkpm_combined_index == NULL) {
 
@@ -7379,7 +7316,6 @@ qmckl_exit_code qmckl_provide_lkpm_combined_index(qmckl_context context)
       ctx->jastrow_champ.lkpm_combined_index = lkpm_combined_index;
     }
 
-    qmckl_exit_code rc;
     rc = qmckl_compute_lkpm_combined_index(context,
                                            ctx->jastrow_champ.cord_num,
                                            ctx->jastrow_champ.dim_c_vector,
@@ -7405,6 +7341,12 @@ qmckl_exit_code qmckl_provide_tmp_c(qmckl_context context)
 
   qmckl_exit_code rc = QMCKL_SUCCESS;
 
+  rc = qmckl_provide_een_rescaled_e(context);
+  if (rc != QMCKL_SUCCESS) return rc;
+
+  rc = qmckl_provide_een_rescaled_n(context);
+  if (rc != QMCKL_SUCCESS) return rc;
+
  /* Compute if necessary */
   if (ctx->date > ctx->jastrow_champ.tmp_c_date) {
 
@@ -7462,6 +7404,12 @@ qmckl_exit_code qmckl_provide_dtmp_c(qmckl_context context)
   qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
   assert (ctx != NULL);
 
+  rc = qmckl_provide_een_rescaled_e_deriv_e(context);
+  if (rc != QMCKL_SUCCESS) return rc;
+
+  rc = qmckl_provide_een_rescaled_n(context);
+  if (rc != QMCKL_SUCCESS) return rc;
+
  /* Compute if necessary */
   if (ctx->date > ctx->jastrow_champ.dtmp_c_date) {
 

org/qmckl_memory.org

@@ -104,12 +104,47 @@ typedef struct qmckl_memory_struct {
   properly freed when the library is de-initialized.
   If the allocation failed, the ~NULL~ pointer is returned.
 
+  The allocated memory block is zeroed using ~memset~.
+
   # Header
   #+begin_src c :tangle (eval h_private_func) :noexport
 void* qmckl_malloc(qmckl_context context,
                    const qmckl_memory_info_struct info);
   #+end_src
 
+  Here's a step-by-step explanation of ~qmckl_malloc~:
+
+  1. The function takes two parameters: a ~qmckl_context~ and a
+ ~qmckl_memory_info_struct~ containing the desired size of the memory
+     block to allocate.
+
+  2. The function checks if the provided ~qmckl_context~ is valid, using the
+ ~qmckl_context_check~ function.
+
+  3. The ~qmckl_context_struct~ pointer is retrieved from the provided
+ ~qmckl_context~.
+
+  4. The function then allocates memory:
+     If the ~HAVE_HPC~ macro is defined, the memory allocation is done using
+     the ~aligned_alloc~ function with a 64-byte alignment, rounding up the
+     requested size to the nearest multiple of 64 bytes. If the ~HAVE_HPC~
+     macro is not defined, the memory allocation is done using the standard
+ ~malloc~ function.
+
+  5  If the allocation fails, the function returns ~NULL~.
+
+  6. The allocated memory block is zeroed using ~memset~.
+
+  7. The function acquires a lock on the ~qmckl_context~ using ~qmckl_lock~.
+
+  8. Inside the locked section, the function checks if the
+ ~qmckl_memory_struct~ is full. If it is, it reallocates a larger array
+     by doubling its size and updating the ~array_size~ member of the
+ ~qmckl_memory_struct~.
+
+  9. The function finds the first available ~qmckl_memory_info_struct~ slot
+     in the element array of the ~qmckl_memory_struct~.
+
   # Source
   #+begin_src c :tangle (eval c)
 void* qmckl_malloc(qmckl_context context, const qmckl_memory_info_struct info) {
@@ -120,6 +155,7 @@ void* qmckl_malloc(qmckl_context context, const qmckl_memory_info_struct info) {
 
  /* Allocate memory and zero it */
 #ifdef HAVE_HPC
+  assert( ((info.size+64) >> 6) << 6 >= info.size );
   void * pointer = aligned_alloc(64, ((info.size+64) >> 6) << 6 );
 #else
   void * pointer = malloc(info.size);
@@ -159,6 +195,7 @@ void* qmckl_malloc(qmckl_context context, const qmckl_memory_info_struct info) {
     memcpy(&(ctx->memory.element[pos]), &info, sizeof(qmckl_memory_info_struct));
     ctx->memory.element[pos].pointer = pointer;
     ctx->memory.n_allocated += (size_t) 1;
+//printf("MALLOC: %5ld  %p\n", ctx->memory.n_allocated, ctx->memory.element[pos].pointer);
   }
   qmckl_unlock(context);
 
@@ -203,6 +240,47 @@ qmckl_exit_code qmckl_free(qmckl_context context,
                            void * const ptr);
   #+end_src
 
+  Here's a step-by-step explanation of the ~qmckl_free~ function:
+
+  1. The function takes two parameters: a ~qmckl_context~ and a pointer to
+     the memory block (~ptr~) that needs to be deallocated.
+
+  2. The function checks if the provided ~qmckl_context~ is valid, using the
+ ~qmckl_context_check~ function. If it is not valid, it returns an error
+     code ~QMCKL_INVALID_CONTEXT~ using the ~qmckl_failwith~ function.
+   
+  3. The function checks if the provided pointer is ~NULL~. If it is, it
+     returns an error code ~QMCKL_INVALID_ARG_2~ using the ~qmckl_failwith~
+     function.
+   
+  4. The ~qmckl_context_struct~ pointer is retrieved from the provided
+ ~qmckl_context~.
+   
+  5. The function acquires a lock on the ~qmckl_context~ using ~qmckl_lock~.
+     
+  6. Inside the locked section, the function searches for the pointer in
+     the element array of the ~qmckl_memory_struct~.
+   
+  7. If the pointer is not found in the array, it releases the lock and
+     returns an error code ~QMCKL_INVALID_ARG_2~ using the ~qmckl_failwith~
+     function.
+     
+  8. If the pointer is found, the memory block is deallocated using the
+     standard ~free~ function.
+
+  9. The ~qmckl_memory_info_struct~ at the found position is zeroed
+     using ~memset~. This marks the slot as available for future
+     allocations.
+     
+ 10. The ~n_allocated~ member of the ~qmckl_memory_struct~ is decremented
+     by one, as the memory block has been deallocated.
+
+ 11. The function releases the lock on the ~qmckl_context~ using ~qmckl_unlock~.
+
+ 12. Finally, the function returns ~QMCKL_SUCCESS~ to indicate
+     successful deallocation of the memory block.
+     
+     
    # Source
    #+begin_src c :tangle (eval c)
 qmckl_exit_code qmckl_free(qmckl_context context, void * const ptr) {
@@ -240,10 +318,13 @@ qmckl_exit_code qmckl_free(qmckl_context context, void * const ptr) {
                             "Pointer not found in context");
     }
 
+ /* Found */
+
     free(ptr);
 
-    memset( &(ctx->memory.element[pos]), 0, sizeof(qmckl_memory_info_struct) );
     ctx->memory.n_allocated -= (size_t) 1;
+//printf("FREE  : %5ld  %p\n", ctx->memory.n_allocated, ctx->memory.element[pos].pointer);
+    ctx->memory.element[pos] = qmckl_memory_info_struct_zero;
   }
   qmckl_unlock(context);