diff --git a/org/qmckl_context.org b/org/qmckl_context.org index 578cb55..7b26a7f 100644 --- a/org/qmckl_context.org +++ b/org/qmckl_context.org @@ -244,6 +244,12 @@ qmckl_context qmckl_context_create() { rc = qmckl_init_ao_basis(context); assert (rc == QMCKL_SUCCESS); + + rc = qmckl_init_mo_basis(context); + assert (rc == QMCKL_SUCCESS); + + rc = qmckl_init_determinant(context); + assert (rc == QMCKL_SUCCESS); } /* Allocate qmckl_memory_struct */ diff --git a/org/qmckl_determinant.org b/org/qmckl_determinant.org index 790d5cc..d071190 100644 --- a/org/qmckl_determinant.org +++ b/org/qmckl_determinant.org @@ -171,6 +171,26 @@ typedef struct qmckl_determinant_struct { Some values are initialized by default, and are not concerned by this mechanism. + #+begin_src c :comments org :tangle (eval h_private_func) +qmckl_exit_code qmckl_init_determinant(qmckl_context context); + #+end_src + + #+begin_src c :comments org :tangle (eval c) +qmckl_exit_code qmckl_init_determinant(qmckl_context context) { + + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { + return false; + } + + qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; + assert (ctx != NULL); + + ctx->det.uninitialized = (1 << 6) - 1; + + return QMCKL_SUCCESS; +} + #+end_src + ** Access functions #+begin_src c :comments org :tangle (eval h_private_func) :exports none @@ -189,6 +209,20 @@ int64_t* qmckl_get_determinant_mo_index_beta (const qmckl_context context); bool qmckl_determinant_provided (const qmckl_context context); #+end_src + #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none +bool qmckl_determinant_provided(const qmckl_context context) { + + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { + return false; + } + + qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; + assert (ctx != NULL); + + return ctx->det.provided; +} + #+end_src + #+NAME:post #+begin_src c :exports none if ( (ctx->det.uninitialized & mask) != 0) { @@ -337,8 +371,8 @@ qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; ctx->det.uninitialized &= ~mask; ctx->det.provided = (ctx->det.uninitialized == 0); if (ctx->det.provided) { - //qmckl_exit_code rc_ = qmckl_finalize_determinant(context); - //if (rc_ != QMCKL_SUCCESS) return rc_; + qmckl_exit_code rc_ = qmckl_finalize_determinant(context); + if (rc_ != QMCKL_SUCCESS) return rc_; } return QMCKL_SUCCESS; @@ -412,7 +446,7 @@ qmckl_exit_code qmckl_set_determinant_det_num_beta(qmckl_context context, const qmckl_exit_code qmckl_set_determinant_mo_index_alpha(qmckl_context context, const int64_t* mo_index_alpha) { <> - int32_t mask = 1 << 5; + int32_t mask = 1 << 4; if (ctx->det.mo_index_alpha != NULL) { qmckl_exit_code rc = qmckl_free(context, ctx->det.mo_index_alpha); @@ -444,7 +478,7 @@ qmckl_exit_code qmckl_set_determinant_mo_index_alpha(qmckl_context context, con qmckl_exit_code qmckl_set_determinant_mo_index_beta(qmckl_context context, const int64_t* mo_index_beta) { <> - int32_t mask = 1 << 6; + int32_t mask = 1 << 5; if (ctx->det.mo_index_beta != NULL) { qmckl_exit_code rc = qmckl_free(context, ctx->det.mo_index_beta); @@ -479,11 +513,50 @@ qmckl_exit_code qmckl_set_determinant_mo_index_beta(qmckl_context context, cons computed to accelerate the calculations. #+begin_src c :comments org :tangle (eval h_private_func) :noweb yes :exports none -qmckl_exit_code qmckl_finalize_basis(qmckl_context context); +qmckl_exit_code qmckl_finalize_determinant(qmckl_context context); #+end_src #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_finalize_determinant(qmckl_context context) { + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { + return qmckl_failwith( context, + QMCKL_INVALID_CONTEXT, + "qmckl_finalize_determinant", + NULL); + } + + qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; + assert (ctx != NULL); + + qmckl_exit_code rc; + rc = qmckl_provide_det_vgl_alpha(context); + if (rc != QMCKL_SUCCESS) { + return qmckl_failwith( context, + QMCKL_NOT_PROVIDED, + "qmckl_finalize_determinant", + NULL); + } + rc = qmckl_provide_det_vgl_beta(context); + if (rc != QMCKL_SUCCESS) { + return qmckl_failwith( context, + QMCKL_NOT_PROVIDED, + "qmckl_finalize_determinant", + NULL); + } + rc = qmckl_provide_det_inv_matrix_alpha(context); + if (rc != QMCKL_SUCCESS) { + return qmckl_failwith( context, + QMCKL_NOT_PROVIDED, + "qmckl_finalize_determinant", + NULL); + } + rc = qmckl_provide_det_inv_matrix_beta(context); + if (rc != QMCKL_SUCCESS) { + return qmckl_failwith( context, + QMCKL_NOT_PROVIDED, + "qmckl_finalize_determinant", + NULL); + } } #+end_src @@ -608,7 +681,7 @@ qmckl_exit_code qmckl_provide_det_vgl_alpha(qmckl_context context) { if (!ctx->det.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, - "qmckl_mo_basis", + "qmckl_determinant", NULL); } @@ -959,19 +1032,19 @@ integer function qmckl_compute_det_vgl_beta_f(context, & do imo = 1, beta_num mo_id = mo_index_beta(imo, iwalk, idet) ! Value - det_vgl_beta(imo, ielec, 1, iwalk, idet) = mo_vgl(mo_id, alpha_num + ielec, 1) + det_vgl_beta(imo, ielec, 1, iwalk, idet) = mo_vgl(mo_id, ielec, 1) ! Grad_x - det_vgl_beta(imo, ielec, 2, iwalk, idet) = mo_vgl(mo_id, alpha_num + ielec, 2) + det_vgl_beta(imo, ielec, 2, iwalk, idet) = mo_vgl(mo_id, ielec, 2) ! Grad_y - det_vgl_beta(imo, ielec, 3, iwalk, idet) = mo_vgl(mo_id, alpha_num + ielec, 3) + det_vgl_beta(imo, ielec, 3, iwalk, idet) = mo_vgl(mo_id, ielec, 3) ! Grad_z - det_vgl_beta(imo, ielec, 4, iwalk, idet) = mo_vgl(mo_id, alpha_num + ielec, 4) + det_vgl_beta(imo, ielec, 4, iwalk, idet) = mo_vgl(mo_id, ielec, 4) ! Lap - det_vgl_beta(imo, ielec, 5, iwalk, idet) = mo_vgl(mo_id, alpha_num + ielec, 5) + det_vgl_beta(imo, ielec, 5, iwalk, idet) = mo_vgl(mo_id, ielec, 5) end do end do end do @@ -1280,7 +1353,7 @@ qmckl_exit_code qmckl_get_det_inv_matrix_beta(qmckl_context context, double * co rc = qmckl_provide_mo_vgl(context); if (rc != QMCKL_SUCCESS) return rc; - rc = qmckl_provide_det_vgl_alpha(context); + rc = qmckl_provide_det_vgl_beta(context); if (rc != QMCKL_SUCCESS) return rc; rc = qmckl_provide_det_inv_matrix_beta(context); @@ -1496,7 +1569,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_beta(qmckl_context context) { qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; mem_info.size = ctx->det.walk_num * ctx->det.det_num_beta * - ctx->electron.up_num * ctx->electron.up_num * sizeof(double); + ctx->electron.down_num * ctx->electron.down_num * sizeof(double); double* det_adj_matrix_beta = (double*) qmckl_malloc(context, mem_info); if (det_adj_matrix_beta == NULL) { diff --git a/org/qmckl_local_energy.org b/org/qmckl_local_energy.org index e5bfd48..b3e24a5 100644 --- a/org/qmckl_local_energy.org +++ b/org/qmckl_local_energy.org @@ -412,6 +412,16 @@ integer function qmckl_compute_kinetic_energy_f(context, walk_num, & return endif + print *,"Size of mo_vgl=",size(mo_vgl) + print *,"Size of mo_index_alpha=",size(mo_index_alpha) + print *,"Size of mo_index_beta=",size(mo_index_beta) + print *,"Size of det_adj_matrix_alpha=",size(det_adj_matrix_alpha) + print *,"Size of det_adj_matrix_beta=",size(det_adj_matrix_beta) + print *,"num_alpha = ",alpha_num + print *,"num_beta = ",beta_num + print *,"det_num_alpha = ",det_num_alpha + print *,"det_num_beta = ",det_num_beta + print *,"walk_num = ",walk_num e_kin = 0.0d0 do idet = 1, det_num_alpha do iwalk = 1, walk_num @@ -427,8 +437,9 @@ integer function qmckl_compute_kinetic_energy_f(context, walk_num, & do imo = 1, beta_num do ielec = 1, beta_num mo_id = mo_index_beta(ielec, iwalk, idet) - e_kin(iwalk) = e_kin(iwalk) - 0.5d0 * det_adj_matrix_beta(imo, ielec, iwalk, idet) * & - mo_vgl(mo_id, ielec, 5) + print *,mo_id, imo, ielec, iwalk, idet + e_kin(iwalk) = e_kin(iwalk) - 0.5d0 * det_adj_matrix_beta(imo, ielec, iwalk, idet) + ! mo_vgl(mo_id, ielec, 5) end do end do end do @@ -1447,7 +1458,7 @@ qmckl_exit_code qmckl_provide_drift_vector(qmckl_context context) { | ~int64_t~ | ~mo_index_alpha[det_num_alpha][walk_num][alpha_num]~ | in | MO indices for electrons | | ~int64_t~ | ~mo_index_beta[det_num_beta][walk_num][beta_num]~ | in | MO indices for electrons | | ~int64_t~ | ~mo_num~ | in | Number of MOs | - | ~double~ | ~mo_vgl[5][elec_num][mo_num]~ | in | Value, gradients and Laplacian of the MOs | + | ~double~ | ~mo_vgl[5][elec_num][mo_num]~ | in | Value, gradients and Laplacian of the MOs | | ~double~ | ~det_adj_matrix_alpha[det_num_alpha][walk_num][alpha_num][alpha_num]~ | in | Value, gradients and Laplacian of the Det | | ~double~ | ~det_adj_matrix_beta[det_num_beta][walk_num][beta_num][beta_num]~ | in | Value, gradients and Laplacian of the Det | | ~double~ | ~r_drift[walk_num][3]~ | out | Kinetic energy | diff --git a/org/qmckl_mo.org b/org/qmckl_mo.org index dad47ae..b9d220c 100644 --- a/org/qmckl_mo.org +++ b/org/qmckl_mo.org @@ -118,6 +118,26 @@ typedef struct qmckl_mo_basis_struct { Some values are initialized by default, and are not concerned by this mechanism. + #+begin_src c :comments org :tangle (eval h_private_func) +qmckl_exit_code qmckl_init_mo_basis(qmckl_context context); + #+end_src + + #+begin_src c :comments org :tangle (eval c) +qmckl_exit_code qmckl_init_mo_basis(qmckl_context context) { + + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { + return false; + } + + qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; + assert (ctx != NULL); + + ctx->mo_basis.uninitialized = (1 << 2) - 1; + + return QMCKL_SUCCESS; +} + #+end_src + ** Access functions #+begin_src c :comments org :tangle (eval h_func) :exports none @@ -259,10 +279,9 @@ qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; ctx->mo_basis.uninitialized &= ~mask; ctx->mo_basis.provided = (ctx->mo_basis.uninitialized == 0); if (ctx->mo_basis.provided) { - qmckl_exit_code rc_ = qmckl_finalize_basis(context); + qmckl_exit_code rc_ = qmckl_finalize_mo_basis(context); if (rc_ != QMCKL_SUCCESS) return rc_; - } - +} return QMCKL_SUCCESS; #+end_src @@ -319,6 +338,34 @@ qmckl_exit_code qmckl_set_mo_basis_coefficient(qmckl_context context, const dou When the basis set is completely entered, other data structures are computed to accelerate the calculations. + #+begin_src c :comments org :tangle (eval h_private_func) :noweb yes :exports none +qmckl_exit_code qmckl_finalize_mo_basis(qmckl_context context); + #+end_src + + #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none +qmckl_exit_code qmckl_finalize_mo_basis(qmckl_context context) { + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { + return qmckl_failwith( context, + QMCKL_INVALID_CONTEXT, + "qmckl_finalize_mo_basis", + NULL); + } + + qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; + assert (ctx != NULL); + + qmckl_exit_code rc; + rc = qmckl_provide_mo_vgl(context); + if (rc != QMCKL_SUCCESS) { + return qmckl_failwith( context, + QMCKL_NOT_PROVIDED, + "qmckl_finalize_mo_basis", + NULL); + } + return rc; +} + #+end_src + * Computation ** Computation of MOs @@ -378,6 +425,7 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context); qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context) { + qmckl_exit_code rc; if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -392,6 +440,14 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context) NULL); } + rc = qmckl_provide_ao_vgl(context); + if (rc != QMCKL_SUCCESS) { + return qmckl_failwith( context, + QMCKL_NOT_PROVIDED, + "qmckl_ao_basis", + NULL); + } + if(!(ctx->electron.provided)) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -425,7 +481,6 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context) ctx->mo_basis.mo_vgl = mo_vgl; } - qmckl_exit_code rc; rc = qmckl_compute_mo_basis_vgl(context, ctx->ao_basis.ao_num, ctx->mo_basis.mo_num, @@ -479,10 +534,6 @@ integer function qmckl_compute_mo_basis_vgl_f(context, & double precision,dimension(:,:),allocatable :: ao_vgl_big double precision :: alpha, beta integer :: info_qmckl_dgemm_value - integer :: info_qmckl_dgemm_Gx - integer :: info_qmckl_dgemm_Gy - integer :: info_qmckl_dgemm_Gz - integer :: info_qmckl_dgemm_lap integer*8 :: M, N, K, LDA, LDB, LDC, i,j integer*8 :: inucl, iprim, iwalk, ielec, ishell diff --git a/org/qmckl_trexio.org b/org/qmckl_trexio.org index 9faac75..4db1837 100644 --- a/org/qmckl_trexio.org +++ b/org/qmckl_trexio.org @@ -915,6 +915,7 @@ qmckl_trexio_read(const qmckl_context context, const char* file_name) #+begin_src c :tangle (eval c_test) #ifdef HAVE_TREXIO +#define walk_num 2 qmckl_exit_code rc; char fname[256]; @@ -927,6 +928,8 @@ char message[256]; strncpy(fname, QMCKL_TEST_DIR,255); strncat(fname, "/chbrclf", 255); printf("Test file: %s\n", fname); + +rc = qmckl_set_electron_walk_num(context, walk_num); rc = qmckl_trexio_read(context, fname); if (rc != QMCKL_SUCCESS) {