From 91946f3ec41d15bb8f637ae63985e44418bd55ed Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Mon, 17 Jan 2022 16:09:41 +0100 Subject: [PATCH] Added size_max to elec_coord --- org/qmckl_ao.org | 6 +- org/qmckl_determinant.org | 134 ++++++++++++++++++------------------- org/qmckl_electron.org | 20 ++++-- org/qmckl_jastrow.org | 42 ++++++------ org/qmckl_local_energy.org | 92 ++++++++++++------------- org/qmckl_mo.org | 40 +++++------ 6 files changed, 173 insertions(+), 161 deletions(-) diff --git a/org/qmckl_ao.org b/org/qmckl_ao.org index 6c04a00..2fd2ed3 100644 --- a/org/qmckl_ao.org +++ b/org/qmckl_ao.org @@ -3312,7 +3312,7 @@ print ( "[7][4][26] : %e"% lf(a,x,y)) assert(qmckl_electron_provided(context)); - rc = qmckl_set_electron_coord (context, 'N', elec_coord); + rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*3*elec_num); assert(rc == QMCKL_SUCCESS); @@ -3726,7 +3726,7 @@ print ( "[1][4][26] : %25.15e"% lf(a,x,y)) assert(qmckl_electron_provided(context)); - rc = qmckl_set_electron_coord (context, 'N', elec_coord); + rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*3*elec_num); assert(rc == QMCKL_SUCCESS); @@ -4752,7 +4752,7 @@ assert (rc == QMCKL_SUCCESS); assert(qmckl_electron_provided(context)); -rc = qmckl_set_electron_coord (context, 'N', elec_coord); +rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*3*elec_num); assert(rc == QMCKL_SUCCESS); diff --git a/org/qmckl_determinant.org b/org/qmckl_determinant.org index 16cafed..e71e2dc 100644 --- a/org/qmckl_determinant.org +++ b/org/qmckl_determinant.org @@ -106,7 +106,7 @@ int main() { | ~mo_index_beta~ | ~mo_index[det_num_beta][walk_num][beta_num]~ | Index of MOs for each walker | Computed data: - + |-----------------------------+------------------------------------------------------+-------------------------------------------------------------------------------------------| | ~up_num~ | ~int64_t~ | Number of number of α electrons | | ~donwn_num~ | ~int64_t~ | Number of number of β electrons | @@ -128,17 +128,17 @@ int main() { | ~det_inv_matrix_beta~ | ~[det_num_beta][walk_num][beta_num][beta_num]~ | Inverse of the β electron slater matrix for each determinant of each walker. | | ~det_inv_matrix_beta_date~ | ~int64_t~ | Date for the Inverse of the β electron slater matrix for each determinant of each walker. | |-----------------------------+------------------------------------------------------+-------------------------------------------------------------------------------------------| - + ** Data structure - + #+begin_src c :comments org :tangle (eval h_private_type) typedef struct qmckl_determinant_struct { char type; - int64_t walk_num; - int64_t det_num_alpha; - int64_t det_num_beta ; - int64_t up_num; - int64_t down_num; + int64_t walk_num; + int64_t det_num_alpha; + int64_t det_num_beta ; + int64_t up_num; + int64_t down_num; int64_t* mo_index_alpha; int64_t* mo_index_beta; @@ -171,10 +171,10 @@ 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) + #+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) { @@ -190,9 +190,9 @@ qmckl_exit_code qmckl_init_determinant(qmckl_context context) { return QMCKL_SUCCESS; } #+end_src - + ** Access functions - + #+begin_src c :comments org :tangle (eval h_private_func) :exports none char qmckl_get_determinant_type (const qmckl_context context); int64_t qmckl_get_determinant_walk_num (const qmckl_context context); @@ -201,7 +201,7 @@ int64_t qmckl_get_determinant_det_num_beta (const qmckl_context context); int64_t* qmckl_get_determinant_mo_index_alpha (const qmckl_context context); int64_t* qmckl_get_determinant_mo_index_beta (const qmckl_context context); #+end_src - + When all the data for the slater determinants have been provided, the following function returns ~true~. @@ -342,7 +342,7 @@ int64_t* qmckl_get_determinant_mo_index_beta (const qmckl_context context) { } #+end_src - + ** Initialization functions To set the basis set, all the following functions need to be @@ -458,7 +458,7 @@ qmckl_exit_code qmckl_set_determinant_mo_index_alpha(qmckl_context context, con } qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; - mem_info.size = ctx->det.walk_num * ctx->det.det_num_alpha * + mem_info.size = ctx->det.walk_num * ctx->det.det_num_alpha * ctx->electron.up_num * sizeof(int64_t); int64_t* new_array = (int64_t*) qmckl_malloc(context, mem_info); if (new_array == NULL) { @@ -490,7 +490,7 @@ qmckl_exit_code qmckl_set_determinant_mo_index_beta(qmckl_context context, cons } qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; - mem_info.size = ctx->det.walk_num * ctx->det.det_num_beta * + mem_info.size = ctx->det.walk_num * ctx->det.det_num_beta * ctx->electron.down_num * sizeof(int64_t); int64_t* new_array = (int64_t*) qmckl_malloc(context, mem_info); if (new_array == NULL) { @@ -572,7 +572,7 @@ qmckl_exit_code qmckl_finalize_determinant(qmckl_context context) { :END: *** Get - + #+begin_src c :comments org :tangle (eval h_func) :noweb yes qmckl_exit_code qmckl_get_det_vgl_alpha(qmckl_context context, double* const det_vgl_alpha); qmckl_exit_code qmckl_get_det_vgl_beta(qmckl_context context, double* const det_vgl_beta); @@ -580,7 +580,7 @@ qmckl_exit_code qmckl_get_det_vgl_beta(qmckl_context context, double* const det_ #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_get_det_vgl_alpha(qmckl_context context, double * const det_vgl_alpha) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -599,7 +599,7 @@ qmckl_exit_code qmckl_get_det_vgl_alpha(qmckl_context context, double * const de qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; assert (ctx != NULL); - size_t sze = 5 * ctx->det.det_num_alpha * ctx->det.walk_num * + size_t sze = 5 * ctx->det.det_num_alpha * ctx->det.walk_num * ctx->electron.up_num * ctx->electron.up_num * sizeof(double); memcpy(det_vgl_alpha, ctx->det.det_vgl_alpha, sze); @@ -607,7 +607,7 @@ qmckl_exit_code qmckl_get_det_vgl_alpha(qmckl_context context, double * const de } qmckl_exit_code qmckl_get_det_vgl_beta(qmckl_context context, double * const det_vgl_beta) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -626,7 +626,7 @@ qmckl_exit_code qmckl_get_det_vgl_beta(qmckl_context context, double * const det qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; assert (ctx != NULL); - size_t sze = 5 * ctx->det.det_num_beta * ctx->det.walk_num * + size_t sze = 5 * ctx->det.det_num_beta * ctx->det.walk_num * ctx->electron.down_num * ctx->electron.down_num * sizeof(double); memcpy(det_vgl_beta, ctx->det.det_vgl_beta, sze); @@ -644,7 +644,7 @@ qmckl_exit_code qmckl_provide_det_vgl_beta(qmckl_context context); #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_provide_det_vgl_alpha(qmckl_context context) { - qmckl_exit_code rc; + qmckl_exit_code rc; if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -658,14 +658,14 @@ qmckl_exit_code qmckl_provide_det_vgl_alpha(qmckl_context context) { "qmckl_electron", NULL); } - + if(!(ctx->electron.provided)) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, "qmckl_electron", NULL); } - + if (!ctx->ao_basis.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -731,7 +731,7 @@ qmckl_exit_code qmckl_provide_det_vgl_alpha(qmckl_context context) { QMCKL_FAILURE, "compute_det_vgl_alpha", "Not yet implemented"); - } + } if (rc != QMCKL_SUCCESS) { return rc; } @@ -757,14 +757,14 @@ qmckl_exit_code qmckl_provide_det_vgl_beta(qmckl_context context) { "qmckl_electron", NULL); } - + if(!(ctx->electron.provided)) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, "qmckl_electron", NULL); } - + if (!ctx->ao_basis.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -806,7 +806,7 @@ qmckl_exit_code qmckl_provide_det_vgl_beta(qmckl_context context) { ctx->det.det_vgl_beta = det_vgl_beta; } - qmckl_exit_code rc; + qmckl_exit_code rc; if (ctx->det.type == 'G') { rc = qmckl_compute_det_vgl_beta(context, ctx->det.det_num_beta, @@ -823,7 +823,7 @@ qmckl_exit_code qmckl_provide_det_vgl_beta(qmckl_context context) { QMCKL_FAILURE, "compute_det_vgl_beta", "Not yet implemented"); - } + } if (rc != QMCKL_SUCCESS) { return rc; } @@ -854,7 +854,7 @@ qmckl_exit_code qmckl_provide_det_vgl_beta(qmckl_context context) { | ~mo_num~ | ~int64_t~ | in | Number of MOs | | ~mo_vgl~ | ~double[5][elec_num][mo_num]~ | in | Value, gradients and Laplacian of the MOs | | ~det_vgl_alpha~ | ~double[det_num_alpha][walk_num][5][alpha_num][alpha_num]~ | out | Value, gradients and Laplacian of the Det | - + #+begin_src f90 :comments org :tangle (eval f) :noweb yes integer function qmckl_compute_det_vgl_alpha_f(context, & det_num_alpha, walk_num, alpha_num, beta_num, elec_num, & @@ -919,7 +919,7 @@ end function qmckl_compute_det_vgl_alpha_f #+end_src #+CALL: generate_c_header(table=qmckl_compute_det_vgl_alpha_args,rettyp=get_value("CRetType"),fname="qmckl_compute_det_vgl_alpha")) - + #+RESULTS: #+begin_src c :tangle (eval h_func) :comments org qmckl_exit_code qmckl_compute_det_vgl_alpha ( @@ -932,7 +932,7 @@ end function qmckl_compute_det_vgl_alpha_f const int64_t* mo_index_alpha, const int64_t mo_num, const double* mo_vgl, - double* const det_vgl_alpha ); + double* const det_vgl_alpha ); #+end_src #+CALL: generate_c_interface(table=qmckl_compute_det_vgl_alpha_args,rettyp=get_value("CRetType"),fname="qmckl_compute_det_vgl_alpha")) @@ -981,7 +981,7 @@ end function qmckl_compute_det_vgl_alpha_f end function qmckl_compute_det_vgl_alpha #+end_src - + *** Compute beta :PROPERTIES: :Name: qmckl_compute_det_vgl_beta @@ -1002,7 +1002,7 @@ end function qmckl_compute_det_vgl_alpha_f | ~mo_num~ | ~int64_t~ | in | Number of MOs | | ~mo_vgl~ | ~double[5][elec_num][mo_num]~ | in | Value, gradients and Laplacian of the MOs | | ~det_vgl_beta~ | ~double[det_num_beta][walk_num][5][beta_num][beta_num]~ | out | Value, gradients and Laplacian of the Det | - + #+begin_src f90 :comments org :tangle (eval f) :noweb yes integer function qmckl_compute_det_vgl_beta_f(context, & det_num_beta, walk_num, alpha_num, beta_num, elec_num, & @@ -1080,9 +1080,9 @@ end function qmckl_compute_det_vgl_beta_f const int64_t* mo_index_beta, const int64_t mo_num, const double* mo_vgl, - double* const det_vgl_beta ); + double* const det_vgl_beta ); #+end_src - + #+CALL: generate_c_interface(table=qmckl_compute_det_vgl_beta_args,rettyp=get_value("CRetType"),fname="qmckl_compute_det_vgl_beta")) #+RESULTS: @@ -1129,7 +1129,7 @@ end function qmckl_compute_det_vgl_beta_f end function qmckl_compute_det_vgl_beta #+end_src - + *** Test #+begin_src c :tangle (eval c_test) :exports none @@ -1146,7 +1146,7 @@ const int64_t nucl_num = chbrclf_nucl_num; const double* nucl_charge = chbrclf_charge; const double* nucl_coord = &(chbrclf_nucl_coord[0][0]); -rc = qmckl_set_electron_num (context, elec_up_num, elec_dn_num); +rc = qmckl_set_electron_num (context, elec_up_num, elec_dn_num); assert (rc == QMCKL_SUCCESS); rc = qmckl_set_electron_walk_num (context, walk_num); @@ -1154,7 +1154,7 @@ assert (rc == QMCKL_SUCCESS); assert(qmckl_electron_provided(context)); -rc = qmckl_set_electron_coord (context, 'N', elec_coord); +rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*elec_num*3); assert(rc == QMCKL_SUCCESS); rc = qmckl_set_nucleus_num (context, nucl_num); @@ -1307,7 +1307,7 @@ rc = qmckl_get_det_vgl_beta(context, &(det_vgl_beta[0][0][0][0][0])); assert (rc == QMCKL_SUCCESS); #+end_src - + ** Inverse of Determinant matrix :PROPERTIES: :Name: qmckl_compute_det_inv_matrix @@ -1316,7 +1316,7 @@ assert (rc == QMCKL_SUCCESS); :END: *** Get - + #+begin_src c :comments org :tangle (eval h_func) :noweb yes qmckl_exit_code qmckl_get_det_inv_matrix_alpha(qmckl_context context, double* const det_inv_matrix_alpha); qmckl_exit_code qmckl_get_det_inv_matrix_beta(qmckl_context context, double* const det_inv_matrix_beta); @@ -1328,7 +1328,7 @@ qmckl_exit_code qmckl_get_det_beta(qmckl_context context, double* const det_adj_ #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_get_det_inv_matrix_alpha(qmckl_context context, double * const det_inv_matrix_alpha) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -1357,7 +1357,7 @@ qmckl_exit_code qmckl_get_det_inv_matrix_alpha(qmckl_context context, double * c } qmckl_exit_code qmckl_get_det_inv_matrix_beta(qmckl_context context, double * const det_inv_matrix_beta) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -1386,7 +1386,7 @@ qmckl_exit_code qmckl_get_det_inv_matrix_beta(qmckl_context context, double * co } qmckl_exit_code qmckl_get_det_adj_matrix_alpha(qmckl_context context, double * const det_adj_matrix_alpha) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -1415,7 +1415,7 @@ qmckl_exit_code qmckl_get_det_adj_matrix_alpha(qmckl_context context, double * c } qmckl_exit_code qmckl_get_det_adj_matrix_beta(qmckl_context context, double * const det_adj_matrix_beta) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -1444,7 +1444,7 @@ qmckl_exit_code qmckl_get_det_adj_matrix_beta(qmckl_context context, double * co } qmckl_exit_code qmckl_get_det_alpha(qmckl_context context, double * const det_value_alpha) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -1473,7 +1473,7 @@ qmckl_exit_code qmckl_get_det_alpha(qmckl_context context, double * const det_va } qmckl_exit_code qmckl_get_det_beta(qmckl_context context, double * const det_value_beta) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -1502,7 +1502,7 @@ qmckl_exit_code qmckl_get_det_beta(qmckl_context context, double * const det_val } #+end_src - + *** Provide #+begin_src c :comments org :tangle (eval h_private_func) :noweb yes :exports none @@ -1526,14 +1526,14 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_alpha(qmckl_context context) { "qmckl_electron", NULL); } - + if(!(ctx->electron.provided)) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, "qmckl_electron", NULL); } - + if (!ctx->ao_basis.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -1562,7 +1562,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_alpha(qmckl_context context) { if (ctx->det.det_inv_matrix_alpha == NULL) { qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; - mem_info.size = ctx->det.walk_num * ctx->det.det_num_alpha * + mem_info.size = ctx->det.walk_num * ctx->det.det_num_alpha * ctx->electron.up_num * ctx->electron.up_num * sizeof(double); double* det_inv_matrix_alpha = (double*) qmckl_malloc(context, mem_info); @@ -1578,7 +1578,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_alpha(qmckl_context context) { if (ctx->det.det_adj_matrix_alpha == NULL) { qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; - mem_info.size = ctx->det.walk_num * ctx->det.det_num_alpha * + mem_info.size = ctx->det.walk_num * ctx->det.det_num_alpha * ctx->electron.up_num * ctx->electron.up_num * sizeof(double); double* det_adj_matrix_alpha = (double*) qmckl_malloc(context, mem_info); @@ -1606,7 +1606,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_alpha(qmckl_context context) { ctx->det.det_value_alpha = det_value_alpha; } - qmckl_exit_code rc; + qmckl_exit_code rc; if (ctx->det.type == 'G') { rc = qmckl_compute_det_inv_matrix_alpha(context, ctx->det.det_num_alpha, @@ -1621,7 +1621,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_alpha(qmckl_context context) { QMCKL_FAILURE, "compute_det_inv_matrix_alpha", "Not yet implemented"); - } + } if (rc != QMCKL_SUCCESS) { return rc; } @@ -1649,14 +1649,14 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_beta(qmckl_context context) { "qmckl_electron", NULL); } - + if(!(ctx->electron.provided)) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, "qmckl_electron", NULL); } - + if (!ctx->ao_basis.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -1685,7 +1685,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_beta(qmckl_context context) { if (ctx->det.det_inv_matrix_beta == NULL) { qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; - mem_info.size = ctx->det.walk_num * ctx->det.det_num_beta * + mem_info.size = ctx->det.walk_num * ctx->det.det_num_beta * ctx->electron.down_num * ctx->electron.down_num * sizeof(double); double* det_inv_matrix_beta = (double*) qmckl_malloc(context, mem_info); @@ -1701,7 +1701,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_beta(qmckl_context context) { if (ctx->det.det_adj_matrix_beta == NULL) { qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; - mem_info.size = ctx->det.walk_num * ctx->det.det_num_beta * + mem_info.size = ctx->det.walk_num * ctx->det.det_num_beta * ctx->electron.down_num * ctx->electron.down_num * sizeof(double); double* det_adj_matrix_beta = (double*) qmckl_malloc(context, mem_info); @@ -1729,7 +1729,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_beta(qmckl_context context) { ctx->det.det_value_beta = det_value_beta; } - qmckl_exit_code rc; + qmckl_exit_code rc; if (ctx->det.type == 'G') { rc = qmckl_compute_det_inv_matrix_beta(context, ctx->det.det_num_beta, @@ -1744,7 +1744,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_beta(qmckl_context context) { QMCKL_FAILURE, "compute_det_inv_matrix_beta", "Not yet implemented"); - } + } if (rc != QMCKL_SUCCESS) { return rc; } @@ -1757,7 +1757,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_beta(qmckl_context context) { return QMCKL_SUCCESS; } #+end_src - + *** Compute alpha :PROPERTIES: :Name: qmckl_compute_det_inv_matrix_alpha @@ -1776,7 +1776,7 @@ qmckl_exit_code qmckl_provide_det_inv_matrix_beta(qmckl_context context) { | ~det_value_alpha~ | ~double[det_num_alpha][walk_num]~ | out | value of determinant matrix | | ~det_adj_matrix_alpha~ | ~double[det_num_alpha][walk_num][alpha_num][alpha_num]~ | out | adjoint of determinant matrix | | ~det_inv_matrix_alpha~ | ~double[det_num_alpha][walk_num][alpha_num][alpha_num]~ | out | inverse of determinant matrix | - + #+begin_src f90 :comments org :tangle (eval f) :noweb yes integer function qmckl_compute_det_inv_matrix_alpha_f(context, & det_num_alpha, walk_num, alpha_num, det_vgl_alpha, det_value_alpha, det_adj_matrix_alpha, det_inv_matrix_alpha) & @@ -1834,7 +1834,7 @@ integer function qmckl_compute_det_inv_matrix_alpha_f(context, & det_inv_matrix_alpha(1:alpha_num, 1:alpha_num, iwalk, idet) = & (1.d0/det_l) * & - det_adj_matrix_alpha(1:alpha_num, 1:alpha_num, iwalk, idet) + det_adj_matrix_alpha(1:alpha_num, 1:alpha_num, iwalk, idet) det_value_alpha(iwalk, idet) = det_l end do @@ -1856,7 +1856,7 @@ end function qmckl_compute_det_inv_matrix_alpha_f const double* det_vgl_alpha, double* const det_value_alpha, double* const det_adj_matrix_alpha, - double* const det_inv_matrix_alpha ); + double* const det_inv_matrix_alpha ); #+end_src #+CALL: generate_c_interface(table=qmckl_det_inv_matrix_alpha_args,rettyp=get_value("CRetType"),fname="qmckl_compute_det_inv_matrix_alpha")) @@ -1918,7 +1918,7 @@ end function qmckl_compute_det_inv_matrix_alpha_f | ~det_value_beta~ | ~double[det_num_beta][walk_num]~ | out | value of determinant matrix | | ~det_adj_matrix_beta~ | ~double[det_num_beta][walk_num][beta_num][beta_num]~ | out | adjoint of determinant matrix | | ~det_inv_matrix_beta~ | ~double[det_num_beta][walk_num][beta_num][beta_num]~ | out | inverse of determinant matrix | - + #+begin_src f90 :comments org :tangle (eval f) :noweb yes integer function qmckl_compute_det_inv_matrix_beta_f(context, & det_num_beta, walk_num, beta_num, det_vgl_beta, det_value_beta, det_adj_matrix_beta, det_inv_matrix_beta) & @@ -1976,7 +1976,7 @@ integer function qmckl_compute_det_inv_matrix_beta_f(context, & det_inv_matrix_beta(1:beta_num, 1:beta_num, iwalk, idet) = & (1.d0/det_l) * & - det_adj_matrix_beta(1:beta_num, 1:beta_num, iwalk, idet) + det_adj_matrix_beta(1:beta_num, 1:beta_num, iwalk, idet) det_value_beta(iwalk, idet) = det_l end do @@ -1999,7 +1999,7 @@ end function qmckl_compute_det_inv_matrix_beta_f const double* det_vgl_beta, double* const det_value_beta, double* const det_adj_matrix_beta, - double* const det_inv_matrix_beta ); + double* const det_inv_matrix_beta ); #+end_src #+CALL: generate_c_interface(table=qmckl_det_inv_matrix_beta_args,rettyp=get_value("CRetType"),fname="qmckl_compute_det_inv_matrix_beta")) diff --git a/org/qmckl_electron.org b/org/qmckl_electron.org index 1f18dcc..89cf86f 100644 --- a/org/qmckl_electron.org +++ b/org/qmckl_electron.org @@ -479,7 +479,7 @@ qmckl_get_electron_coord (const qmckl_context context, const char transp, double #+begin_src c :comments org :tangle (eval h_func) qmckl_exit_code qmckl_set_electron_num (qmckl_context context, const int64_t up_num, const int64_t down_num); qmckl_exit_code qmckl_set_electron_walk_num (qmckl_context context, const int64_t walk_num); -qmckl_exit_code qmckl_set_electron_coord (qmckl_context context, const char transp, const double* coord); +qmckl_exit_code qmckl_set_electron_coord (qmckl_context context, const char transp, const double* coord, const int64_t size_max); qmckl_exit_code qmckl_set_electron_rescale_factor_ee (qmckl_context context, const double kappa_ee); qmckl_exit_code qmckl_set_electron_rescale_factor_en (qmckl_context context, const double kappa_en); @@ -664,7 +664,11 @@ end interface #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code -qmckl_set_electron_coord(qmckl_context context, const char transp, const double* coord) { +qmckl_set_electron_coord(qmckl_context context, + const char transp, + const double* coord, + const int64_t size_max) +{ <> @@ -705,6 +709,13 @@ qmckl_set_electron_coord(qmckl_context context, const char transp, const double* "walk_num is not set"); } + if (size_max < walk_num*3*elec_num) { + return qmckl_failwith( context, + QMCKL_INVALID_ARG_4, + "qmckl_set_electron_coord", + "destination array is too small"); + } + /* If num and walk_num are set, the arrays should be allocated */ assert (ctx->electron.coord_old != NULL); assert (ctx->electron.coord_new != NULL); @@ -742,7 +753,7 @@ qmckl_set_electron_coord(qmckl_context context, const char transp, const double* #+begin_src f90 :comments org :tangle (eval fh_func) :noweb yes interface - integer(c_int32_t) function qmckl_set_electron_coord(context, transp, coord) bind(C) + integer(c_int32_t) function qmckl_set_electron_coord(context, transp, coord, size_max) bind(C) use, intrinsic :: iso_c_binding import implicit none @@ -750,6 +761,7 @@ interface integer (c_int64_t) , intent(in) , value :: context character , intent(in) , value :: transp double precision , intent(in) :: coord(*) + integer (c_int64_t) , intent(in) , value :: size_max end function end interface #+end_src @@ -848,7 +860,7 @@ assert(w == walk_num); assert(qmckl_electron_provided(context)); -rc = qmckl_set_electron_coord (context, 'N', elec_coord); +rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*elec_num*3); assert(rc == QMCKL_SUCCESS); double elec_coord2[walk_num*3*elec_num]; diff --git a/org/qmckl_jastrow.org b/org/qmckl_jastrow.org index 0c43700..dd0eea3 100644 --- a/org/qmckl_jastrow.org +++ b/org/qmckl_jastrow.org @@ -103,26 +103,26 @@ int main() { computed data: - | Variable | Type | In/Out | Description | - |------------+-----------------------------------------------------------------------+---------------------------------------------------------------------------------------------------------+-------------| - | ~dim_cord_vect~ | ~int64_t~ | Number of unique C coefficients | | - | ~dim_cord_vect_date~ | ~uint64_t~ | Number of unique C coefficients | | - | ~asymp_jasb~ | ~double[2]~ | Asymptotic component | | - | ~asymp_jasb_date~ | ~uint64_t~ | Asymptotic component | | - | ~cord_vect_full~ | ~double[dim_cord_vect][nucl_num]~ | vector of non-zero coefficients | | - | ~cord_vect_full_date~ | ~uint64_t~ | Keep track of changes here | | - | ~lkpm_combined_index~ | ~int64_t[4][dim_cord_vect]~ | Transform l,k,p, and m into consecutive indices | | - | ~lkpm_combined_index_date~ | ~uint64_t~ | Transform l,k,p, and m into consecutive indices | | - | ~tmp_c~ | ~double[elec_num][nucl_num][ncord + 1][ncord][walk_num]~ | vector of non-zero coefficients | | - | ~dtmp_c~ | ~double[elec_num][4][nucl_num][ncord + 1][ncord][walk_num]~ | vector of non-zero coefficients | | - | ~een_rescaled_e~ | ~double[walk_num][elec_num][elec_num][0:cord_num]~ | The electron-electron rescaled distances raised to the powers defined by cord | | - | ~een_rescaled_e_date~ | ~uint64_t~ | Keep track of the date of creation | | - | ~een_rescaled_n~ | ~double[walk_num][elec_num][nucl_num][0:cord_num]~ | The electron-electron rescaled distances raised to the powers defined by cord | | - | ~een_rescaled_n_date~ | ~uint64_t~ | Keep track of the date of creation | | - | ~een_rescaled_e_deriv_e~ | ~double[walk_num][elec_num][4][elec_num][0:cord_num]~ | The electron-electron rescaled distances raised to the powers defined by cord derivatives wrt electrons | | - | ~een_rescaled_e_deriv_e_date~ | ~uint64_t~ | Keep track of the date of creation | | - | ~een_rescaled_n_deriv_e~ | ~double[walk_num][elec_num][4][nucl_num][0:cord_num]~ | The electron-electron rescaled distances raised to the powers defined by cord derivatives wrt electrons | | - | ~een_rescaled_n_deriv_e_date~ | ~uint64_t~ | Keep track of the date of creation | | + | Variable | Type | In/Out | Description | + |-------------------------------+-------------------------------------------------------------+---------------------------------------------------------------------------------------------------------+-------------| + | ~dim_cord_vect~ | ~int64_t~ | Number of unique C coefficients | | + | ~dim_cord_vect_date~ | ~uint64_t~ | Number of unique C coefficients | | + | ~asymp_jasb~ | ~double[2]~ | Asymptotic component | | + | ~asymp_jasb_date~ | ~uint64_t~ | Asymptotic component | | + | ~cord_vect_full~ | ~double[dim_cord_vect][nucl_num]~ | vector of non-zero coefficients | | + | ~cord_vect_full_date~ | ~uint64_t~ | Keep track of changes here | | + | ~lkpm_combined_index~ | ~int64_t[4][dim_cord_vect]~ | Transform l,k,p, and m into consecutive indices | | + | ~lkpm_combined_index_date~ | ~uint64_t~ | Transform l,k,p, and m into consecutive indices | | + | ~tmp_c~ | ~double[elec_num][nucl_num][ncord + 1][ncord][walk_num]~ | vector of non-zero coefficients | | + | ~dtmp_c~ | ~double[elec_num][4][nucl_num][ncord + 1][ncord][walk_num]~ | vector of non-zero coefficients | | + | ~een_rescaled_e~ | ~double[walk_num][elec_num][elec_num][0:cord_num]~ | The electron-electron rescaled distances raised to the powers defined by cord | | + | ~een_rescaled_e_date~ | ~uint64_t~ | Keep track of the date of creation | | + | ~een_rescaled_n~ | ~double[walk_num][elec_num][nucl_num][0:cord_num]~ | The electron-electron rescaled distances raised to the powers defined by cord | | + | ~een_rescaled_n_date~ | ~uint64_t~ | Keep track of the date of creation | | + | ~een_rescaled_e_deriv_e~ | ~double[walk_num][elec_num][4][elec_num][0:cord_num]~ | The electron-electron rescaled distances raised to the powers defined by cord derivatives wrt electrons | | + | ~een_rescaled_e_deriv_e_date~ | ~uint64_t~ | Keep track of the date of creation | | + | ~een_rescaled_n_deriv_e~ | ~double[walk_num][elec_num][4][nucl_num][0:cord_num]~ | The electron-electron rescaled distances raised to the powers defined by cord derivatives wrt electrons | | + | ~een_rescaled_n_deriv_e_date~ | ~uint64_t~ | Keep track of the date of creation | | For H2O we have the following data: @@ -1088,7 +1088,7 @@ assert(w == walk_num); assert(qmckl_electron_provided(context)); -rc = qmckl_set_electron_coord (context, 'N', elec_coord); +rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*3*elec_num); assert(rc == QMCKL_SUCCESS); double elec_coord2[walk_num*3*elec_num]; diff --git a/org/qmckl_local_energy.org b/org/qmckl_local_energy.org index c3c2605..c29974d 100644 --- a/org/qmckl_local_energy.org +++ b/org/qmckl_local_energy.org @@ -13,7 +13,7 @@ E_L = KE + PE Where the kinetic energy is given as: \[ -KE = -\frac{1}{2}\frac{\bigtriangleup \Psi}{\Psi} +KE = -\frac{1}{2}\frac{\bigtriangleup \Psi}{\Psi} \] The laplacian of the wavefunction in the single-determinant @@ -59,7 +59,7 @@ F(x) = 2\frac{\nabla \Psi(r)}{\Psi(r)} y = x + D F(x) \delta t + \chi \] -Where \[\chi\] is a random number with gaussian distribution centered at 0 +Where \[\chi\] is a random number with gaussian distribution centered at 0 and width of \[2D\delta t\]. Here \[D\] is the drift parameter. 3. Acceptance probability - \[min\left[1, q(y,x)\right]\] @@ -102,7 +102,7 @@ int main() { qmckl_exit_code rc; #+end_src - + #+begin_src c :tangle (eval c) #ifdef HAVE_CONFIG_H #include "config.h" @@ -138,7 +138,7 @@ int main() { | | Computed data: - + |--------------+---------------------------+----------------------------| | ~e_kin~ | ~[walk_num]~ | total kinetic energy | | ~e_pot~ | ~[walk_num]~ | total potential energy | @@ -147,9 +147,9 @@ int main() { | ~y_move~ | ~[3][walk_num]~ | The diffusion move | | ~accep_prob~ | ~[walk_num]~ | The acceptance probability | |--------------+---------------------------+----------------------------| - + ** Data structure - + #+begin_src c :comments org :tangle (eval h_private_type) typedef struct qmckl_local_energy_struct { double * e_kin; @@ -188,7 +188,7 @@ typedef struct qmckl_local_energy_struct { Where the kinetic energy is given as: \[ -KE = -\frac{1}{2}\frac{\bigtriangleup \Psi}{\Psi} +KE = -\frac{1}{2}\frac{\bigtriangleup \Psi}{\Psi} \] The laplacian of the wavefunction in the single-determinant @@ -199,14 +199,14 @@ case is given as follows: \] *** Get - + #+begin_src c :comments org :tangle (eval h_func) :noweb yes qmckl_exit_code qmckl_get_kinetic_energy(qmckl_context context, double* const kinetic_energy); #+end_src #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_get_kinetic_energy(qmckl_context context, double * const kinetic_energy) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -245,7 +245,7 @@ qmckl_exit_code qmckl_provide_kinetic_energy(qmckl_context context); #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_provide_kinetic_energy(qmckl_context context) { - qmckl_exit_code rc; + qmckl_exit_code rc; if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -259,14 +259,14 @@ qmckl_exit_code qmckl_provide_kinetic_energy(qmckl_context context) { "qmckl_electron", NULL); } - + if(!(ctx->electron.provided)) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, "qmckl_electron", NULL); } - + if (!ctx->ao_basis.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -344,7 +344,7 @@ qmckl_exit_code qmckl_provide_kinetic_energy(qmckl_context context) { QMCKL_FAILURE, "compute_kinetic_energy", "Not yet implemented"); - } + } if (rc != QMCKL_SUCCESS) { return rc; } @@ -380,7 +380,7 @@ qmckl_exit_code qmckl_provide_kinetic_energy(qmckl_context context) { | ~double~ | ~det_inv_matrix_alpha[det_num_alpha][walk_num][alpha_num][alpha_num]~ | in | Value, gradients and Laplacian of the Det | | ~double~ | ~det_inv_matrix_beta[det_num_beta][walk_num][beta_num][beta_num]~ | in | Value, gradients and Laplacian of the Det | | ~double~ | ~e_kin[walk_num]~ | out | Kinetic energy | - + #+begin_src f90 :comments org :tangle (eval f) :noweb yes integer function qmckl_compute_kinetic_energy_f(context, walk_num, & det_num_alpha, det_num_beta, alpha_num, beta_num, elec_num, mo_index_alpha, mo_index_beta, & @@ -493,7 +493,7 @@ end function qmckl_compute_kinetic_energy_f const double* det_value_beta, const double* det_inv_matrix_alpha, const double* det_inv_matrix_beta, - double* const e_kin ); + double* const e_kin ); #+end_src #+CALL: generate_c_interface(table=qmckl_compute_kinetic_energy_args,rettyp=get_value("CRetType"),fname="qmckl_compute_kinetic_energy")) @@ -560,7 +560,7 @@ end function qmckl_compute_kinetic_energy_f end function qmckl_compute_kinetic_energy #+end_src - + *** Test #+begin_src c :tangle (eval c_test) :exports none @@ -576,7 +576,7 @@ const int64_t nucl_num = chbrclf_nucl_num; const double* nucl_charge = chbrclf_charge; const double* nucl_coord = &(chbrclf_nucl_coord[0][0]); -rc = qmckl_set_electron_num (context, elec_up_num, elec_dn_num); +rc = qmckl_set_electron_num (context, elec_up_num, elec_dn_num); assert (rc == QMCKL_SUCCESS); rc = qmckl_set_electron_walk_num (context, walk_num); @@ -584,7 +584,7 @@ assert (rc == QMCKL_SUCCESS); assert(qmckl_electron_provided(context)); -rc = qmckl_set_electron_coord (context, 'N', elec_coord); +rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*elec_num*3); assert(rc == QMCKL_SUCCESS); rc = qmckl_set_nucleus_num (context, nucl_num); @@ -756,7 +756,7 @@ rc = qmckl_get_kinetic_energy(context, &(kinetic_energy[0])); assert (rc == QMCKL_SUCCESS); #+end_src - + ** Potential energy :PROPERTIES: :Name: qmckl_compute_potential_energy @@ -786,14 +786,14 @@ contributions. \] *** Get - + #+begin_src c :comments org :tangle (eval h_func) :noweb yes qmckl_exit_code qmckl_get_potential_energy(qmckl_context context, double* const potential_energy); #+end_src #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_get_potential_energy(qmckl_context context, double * const potential_energy) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -839,21 +839,21 @@ qmckl_exit_code qmckl_provide_potential_energy(qmckl_context context) { qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; assert (ctx != NULL); - qmckl_exit_code rc; + qmckl_exit_code rc; if(!(ctx->nucleus.provided)) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, "qmckl_electron", NULL); } - + if(!(ctx->electron.provided)) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, "qmckl_electron", NULL); } - + if (!ctx->ao_basis.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -932,7 +932,7 @@ qmckl_exit_code qmckl_provide_potential_energy(qmckl_context context) { QMCKL_FAILURE, "compute_potential_energy", "Not yet implemented"); - } + } if (rc != QMCKL_SUCCESS) { return rc; } @@ -960,7 +960,7 @@ qmckl_exit_code qmckl_provide_potential_energy(qmckl_context context) { | ~double~ | ~en_pot[walk_num]~ | in | en potential | | ~double~ | ~repulsion~ | in | en potential | | ~double~ | ~e_pot[walk_num]~ | out | Potential energy | - + #+begin_src f90 :comments org :tangle (eval f) :noweb yes integer function qmckl_compute_potential_energy_f(context, walk_num, & elec_num, nucl_num, ee_pot, en_pot, repulsion, e_pot) & @@ -1014,7 +1014,7 @@ end function qmckl_compute_potential_energy_f const double* ee_pot, const double* en_pot, const double repulsion, - double* const e_pot ); + double* const e_pot ); #+end_src #+CALL: generate_c_interface(table=qmckl_compute_potential_energy_args,rettyp=get_value("CRetType"),fname="qmckl_compute_potential_energy")) @@ -1043,7 +1043,7 @@ end function qmckl_compute_potential_energy_f end function qmckl_compute_potential_energy #+end_src - + *** Test #+begin_src c :tangle (eval c_test) :exports none // Calculate the Potential energy @@ -1070,14 +1070,14 @@ E_L = KE + PE *** Get - + #+begin_src c :comments org :tangle (eval h_func) :noweb yes qmckl_exit_code qmckl_get_local_energy(qmckl_context context, double* const local_energy); #+end_src #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_get_local_energy(qmckl_context context, double * const local_energy) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -1129,14 +1129,14 @@ qmckl_exit_code qmckl_provide_local_energy(qmckl_context context) { "qmckl_electron", NULL); } - + if(!(ctx->electron.provided)) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, "qmckl_electron", NULL); } - + if (!ctx->ao_basis.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -1205,7 +1205,7 @@ qmckl_exit_code qmckl_provide_local_energy(qmckl_context context) { QMCKL_FAILURE, "compute_local_energy", "Not yet implemented"); - } + } if (rc != QMCKL_SUCCESS) { return rc; } @@ -1230,7 +1230,7 @@ qmckl_exit_code qmckl_provide_local_energy(qmckl_context context) { | ~double~ | ~e_kin[walk_num]~ | in | e kinetic | | ~double~ | ~e_pot[walk_num]~ | in | e potential | | ~double~ | ~e_local[walk_num]~ | out | local energy | - + #+begin_src f90 :comments org :tangle (eval f) :noweb yes integer function qmckl_compute_local_energy_f(context, walk_num, & e_kin, e_pot, e_local) & @@ -1273,7 +1273,7 @@ end function qmckl_compute_local_energy_f const int64_t walk_num, const double* e_kin, const double* e_pot, - double* const e_local ); + double* const e_local ); #+end_src #+CALL: generate_c_interface(table=qmckl_compute_local_energy_args,rettyp=get_value("CRetType"),fname="qmckl_compute_local_energy")) @@ -1299,7 +1299,7 @@ end function qmckl_compute_local_energy_f end function qmckl_compute_local_energy #+end_src - + *** Test #+begin_src c :tangle (eval c_test) :exports none // Calculate the Local energy @@ -1318,7 +1318,7 @@ assert (rc == QMCKL_SUCCESS); :FRetType: qmckl_exit_code :END: -The drift vector is calculated as the ration of the gradient +The drift vector is calculated as the ration of the gradient with the determinant of the wavefunction. \[ @@ -1326,14 +1326,14 @@ with the determinant of the wavefunction. \] *** Get - + #+begin_src c :comments org :tangle (eval h_func) :noweb yes qmckl_exit_code qmckl_get_drift_vector(qmckl_context context, double* const drift_vector); #+end_src #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_get_drift_vector(qmckl_context context, double * const drift_vector) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -1385,14 +1385,14 @@ qmckl_exit_code qmckl_provide_drift_vector(qmckl_context context) { "qmckl_electron", NULL); } - + if(!(ctx->electron.provided)) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, "qmckl_electron", NULL); } - + if (!ctx->ao_basis.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -1433,7 +1433,7 @@ qmckl_exit_code qmckl_provide_drift_vector(qmckl_context context) { ctx->local_energy.r_drift = r_drift; } - qmckl_exit_code rc; + qmckl_exit_code rc; if (ctx->det.type == 'G') { rc = qmckl_compute_drift_vector(context, ctx->det.walk_num, @@ -1454,7 +1454,7 @@ qmckl_exit_code qmckl_provide_drift_vector(qmckl_context context) { QMCKL_FAILURE, "compute_drift_vector", "Not yet implemented"); - } + } if (rc != QMCKL_SUCCESS) { return rc; } @@ -1488,7 +1488,7 @@ qmckl_exit_code qmckl_provide_drift_vector(qmckl_context context) { | ~double~ | ~det_inv_matrix_alpha[det_num_alpha][walk_num][alpha_num][alpha_num]~ | in | Value, gradients and Laplacian of the Det | | ~double~ | ~det_inv_matrix_beta[det_num_beta][walk_num][beta_num][beta_num]~ | in | Value, gradients and Laplacian of the Det | | ~double~ | ~r_drift[walk_num][elec_num][3]~ | out | Kinetic energy | - + #+begin_src f90 :comments org :tangle (eval f) :noweb yes integer function qmckl_compute_drift_vector_f(context, walk_num, & det_num_alpha, det_num_beta, alpha_num, beta_num, elec_num, mo_index_alpha, mo_index_beta, & @@ -1593,7 +1593,7 @@ end function qmckl_compute_drift_vector_f const double* mo_vgl, const double* det_inv_matrix_alpha, const double* det_inv_matrix_beta, - double* const r_drift ); + double* const r_drift ); #+end_src #+CALL: generate_c_interface(table=qmckl_compute_drift_vector_args,rettyp=get_value("CRetType"),fname="qmckl_compute_drift_vector")) @@ -1654,7 +1654,7 @@ end function qmckl_compute_drift_vector_f end function qmckl_compute_drift_vector #+end_src - + *** Test #+begin_src c :tangle (eval c_test) :exports none // Calculate the Drift vector diff --git a/org/qmckl_mo.org b/org/qmckl_mo.org index 1da14b1..a07bba4 100644 --- a/org/qmckl_mo.org +++ b/org/qmckl_mo.org @@ -2,7 +2,7 @@ #+SETUPFILE: ../tools/theme.setup #+INCLUDE: ../tools/lib.org -The molecular orbitals (MOs) are defined in the basis of AOs along with a AO to MO +The molecular orbitals (MOs) are defined in the basis of AOs along with a AO to MO coefficient matrix \[C\]. Using these coefficients (e.g. from Hartree Fock SCF method) the MOs are defined as follows: @@ -89,7 +89,7 @@ int main() { | ~coefficient~ | ~[mo_num][ao_num]~ | Orbital coefficients | Computed data: - + |---------------+-------------------------+----------------------------------------------------------------------------------------| |---------------+-------------------------+----------------------------------------------------------------------------------------| | ~mo_vgl~ | ~[5][elec_num][mo_num]~ | Value, gradients, Laplacian of the MOs at electron positions | @@ -100,7 +100,7 @@ int main() { #+begin_src c :comments org :tangle (eval h_private_type) typedef struct qmckl_mo_basis_struct { - int64_t mo_num; + int64_t mo_num; double * coefficient; double * mo_vgl; @@ -118,10 +118,10 @@ 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) + #+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) { @@ -137,7 +137,7 @@ qmckl_exit_code qmckl_init_mo_basis(qmckl_context context) { return QMCKL_SUCCESS; } #+end_src - + ** Access functions #+begin_src c :comments org :tangle (eval h_func) :exports none @@ -372,7 +372,7 @@ qmckl_exit_code qmckl_get_mo_basis_vgl(qmckl_context context, double* const mo_v #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_get_mo_basis_vgl(qmckl_context context, double* const mo_vgl) { - + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -419,7 +419,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; + qmckl_exit_code rc; if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_NULL_CONTEXT; } @@ -448,7 +448,7 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context) "qmckl_electron", NULL); } - + if (!ctx->mo_basis.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -492,7 +492,7 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context) return QMCKL_SUCCESS; } #+end_src - + *** Compute :PROPERTIES: :Name: qmckl_compute_mo_basis_vgl @@ -509,7 +509,7 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context) | ~double~ | ~ao_vgl[5][elec_num][ao_num]~ | in | Value, gradients and Laplacian of the AOs | | ~double~ | ~mo_vgl[5][elec_num][mo_num]~ | out | Value, gradients and Laplacian of the MOs | - + #+begin_src f90 :comments org :tangle (eval f) :noweb yes integer function qmckl_compute_mo_basis_vgl_f(context, & ao_num, mo_num, elec_num, & @@ -547,7 +547,7 @@ integer function qmckl_compute_mo_basis_vgl_f(context, & info = QMCKL_SUCCESS info_qmckl_dgemm_value = QMCKL_SUCCESS - + ! Don't compute exponentials when the result will be almost zero. ! TODO : Use numerical precision here cutoff = -dlog(1.d-15) @@ -571,7 +571,7 @@ integer function qmckl_compute_mo_basis_vgl_f(context, & !end do ao_vgl_big = reshape(ao_vgl(:, :, :),(/ao_num, elec_num*5_8/)) - LDB = size(ao_vgl_big,1) + LDB = size(ao_vgl_big,1) LDC = size(mo_vgl_big,1) info = qmckl_dgemm(context,TransA, TransB, M, N, K, alpha, & @@ -598,7 +598,7 @@ end function qmckl_compute_mo_basis_vgl_f const int64_t elec_num, const double* coef_normalized, const double* ao_vgl, - double* const mo_vgl ); + double* const mo_vgl ); #+end_src @@ -634,7 +634,7 @@ end function qmckl_compute_mo_basis_vgl_f import numpy as np def f(a,x,y): - return np.sum( [c * np.exp( -b*(np.linalg.norm(x-y))**2) for b,c in a] ) + return np.sum( [c * np.exp( -b*(np.linalg.norm(x-y))**2) for b,c in a] ) def df(a,x,y,n): h0 = 1.e-6 @@ -706,7 +706,7 @@ const int64_t nucl_num = chbrclf_nucl_num; const double* nucl_charge = chbrclf_charge; const double* nucl_coord = &(chbrclf_nucl_coord[0][0]); -rc = qmckl_set_electron_num (context, elec_up_num, elec_dn_num); +rc = qmckl_set_electron_num (context, elec_up_num, elec_dn_num); assert (rc == QMCKL_SUCCESS); rc = qmckl_set_electron_walk_num (context, walk_num); @@ -714,7 +714,7 @@ assert (rc == QMCKL_SUCCESS); assert(qmckl_electron_provided(context)); -rc = qmckl_set_electron_coord (context, 'N', elec_coord); +rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*elec_num*3); assert(rc == QMCKL_SUCCESS); rc = qmckl_set_nucleus_num (context, nucl_num); @@ -848,7 +848,7 @@ assert (rc == QMCKL_SUCCESS); // elec_coord[0] = point_x[i]; // elec_coord[1] = point_y[j]; // elec_coord[2] = point_z[k]; -// rc = qmckl_set_electron_coord (context, 'N', elec_coord); +// rc = qmckl_set_electron_coord (context, 'N', elec_coord); // assert(rc == QMCKL_SUCCESS); // // // Calculate value of MO (1st electron) @@ -876,8 +876,8 @@ printf(" mo_vgl mo_vgl[1][26][223] %25.15e\n", mo_vgl[1][2][3]); printf(" mo_vgl mo_vgl[0][26][224] %25.15e\n", mo_vgl[0][2][3]); printf(" mo_vgl mo_vgl[1][26][224] %25.15e\n", mo_vgl[1][2][3]); printf("\n"); -} - +} + #+end_src * End of files :noexport: