diff --git a/org/qmckl_electron.org b/org/qmckl_electron.org index 284f3ef..8e9a298 100644 --- a/org/qmckl_electron.org +++ b/org/qmckl_electron.org @@ -63,29 +63,29 @@ int main() { The following data stored in the context: - | ~uninitialized~ | ~int32_t~ | Keeps bit set for uninitialized data | - | ~num~ | ~int64_t~ | Total number of electrons | - | ~up_num~ | ~int64_t~ | Number of up-spin electrons | - | ~down_num~ | ~int64_t~ | Number of down-spin electrons | - | ~walk_num~ | ~int64_t~ | Number of walkers | - | ~rescale_factor_kappa_ee~ | ~double~ | The distance scaling factor | - | ~rescale_factor_kappa_en~ | ~double~ | The distance scaling factor | - | ~provided~ | ~bool~ | If true, ~electron~ is valid | - | ~coord_new~ | ~double[walk_num][3][num]~ | New set of electron coordinates | - | ~coord_old~ | ~double[walk_num][3][num]~ | Old set of electron coordinates | - | ~coord_new_date~ | ~uint64_t~ | Last modification date of the coordinates | - | ~ee_distance~ | ~double[walk_num][num][num]~ | Electron-electron distances | - | ~ee_distance_date~ | ~uint64_t~ | Last modification date of the electron-electron distances | - | ~en_distance~ | ~double[walk_num][nucl_num][num]~ | Electron-nucleus distances | - | ~en_distance_date~ | ~uint64_t~ | Last modification date of the electron-electron distances | - | ~ee_distance_rescaled~ | ~double[walk_num][num][num]~ | Electron-electron rescaled distances | - | ~ee_distance_rescaled_date~ | ~uint64_t~ | Last modification date of the electron-electron distances | - | ~ee_distance_rescaled_deriv_e~ | ~double[walk_num][4][num][num]~ | Electron-electron rescaled distances derivatives | - | ~ee_distance_rescaled_deriv_e_date~ | ~uint64_t~ | Last modification date of the electron-electron distance derivatives | - | ~en_distance_rescaled~ | ~double[walk_num][nucl_num][num]~ | Electron-nucleus distances | - | ~en_distance_rescaled_date~ | ~uint64_t~ | Last modification date of the electron-electron distances | - | ~en_distance_rescaled_deriv_e~ | ~double[walk_num][4][num][num]~ | Electron-electron rescaled distances derivatives | - | ~en_distance_rescaled_deriv_e_date~ | ~uint64_t~ | Last modification date of the electron-electron distance derivatives | + | ~uninitialized~ | ~int32_t~ | Keeps bit set for uninitialized data | + | ~num~ | ~int64_t~ | Total number of electrons | + | ~up_num~ | ~int64_t~ | Number of up-spin electrons | + | ~down_num~ | ~int64_t~ | Number of down-spin electrons | + | ~walk_num~ | ~int64_t~ | Number of walkers | + | ~rescale_factor_kappa_ee~ | ~double~ | The distance scaling factor | + | ~rescale_factor_kappa_en~ | ~double~ | The distance scaling factor | + | ~provided~ | ~bool~ | If true, ~electron~ is valid | + | ~coord_new~ | ~double[walk_num][3][num]~ | New set of electron coordinates | + | ~coord_old~ | ~double[walk_num][3][num]~ | Old set of electron coordinates | + | ~coord_new_date~ | ~uint64_t~ | Last modification date of the coordinates | + | ~ee_distance~ | ~double[walk_num][num][num]~ | Electron-electron distances | + | ~ee_distance_date~ | ~uint64_t~ | Last modification date of the electron-electron distances | + | ~en_distance~ | ~double[walk_num][nucl_num][num]~ | Electron-nucleus distances | + | ~en_distance_date~ | ~uint64_t~ | Last modification date of the electron-electron distances | + | ~ee_distance_rescaled~ | ~double[walk_num][num][num]~ | Electron-electron rescaled distances | + | ~ee_distance_rescaled_date~ | ~uint64_t~ | Last modification date of the electron-electron distances | + | ~ee_distance_rescaled_deriv_e~ | ~double[walk_num][4][num][num]~ | Electron-electron rescaled distances derivatives | + | ~ee_distance_rescaled_deriv_e_date~ | ~uint64_t~ | Last modification date of the electron-electron distance derivatives | + | ~en_distance_rescaled~ | ~double[walk_num][nucl_num][num]~ | Electron-nucleus distances | + | ~en_distance_rescaled_date~ | ~uint64_t~ | Last modification date of the electron-electron distances | + | ~en_distance_rescaled_deriv_e~ | ~double[walk_num][4][nucl_num][num]~ | Electron-electron rescaled distances derivatives | + | ~en_distance_rescaled_deriv_e_date~ | ~uint64_t~ | Last modification date of the electron-electron distance derivatives | ** Data structure @@ -2112,7 +2112,7 @@ assert(fabs(en_distance_rescaled[1][0][1] - 0.9584331688679852) < 1.e-12); The rescaled distances which is given as $R = (1 - \exp{-\kappa r})/\kappa$ needs to be perturbed with respect to the nuclear coordinates. This data is stored in the ~en_distance_rescaled_deriv_e~ tensor. The - The first three elements of this three index tensor ~[4][num][num]~ gives the + The first three elements of this three index tensor ~[4][nucl_num][elec_num]~ gives the derivatives in the x, y, and z directions $dx, dy, dz$ and the last index gives the Laplacian $\partial x^2 + \partial y^2 + \partial z^2$. diff --git a/org/qmckl_jastrow.org b/org/qmckl_jastrow.org index 461a445..d0da070 100644 --- a/org/qmckl_jastrow.org +++ b/org/qmckl_jastrow.org @@ -2129,7 +2129,6 @@ assert(fabs(factor_ee_deriv_e[0][3][0]-1.5111672803213185 ) < 1.e-12); #+end_src - ** Electron-nucleus component \(f_{en}\) Calculate the electron-electron jastrow component ~factor_en~ using the ~aord_vector~ @@ -2422,6 +2421,333 @@ assert(fabs(factor_en[0]+5.865822569188727) < 1.e-12); #+end_src +** Electron-nucleus component derivative \(f'_{en}\) + + Calculate the electron-electron jastrow component ~factor_en_deriv_e~ derivative + with respect to the electron coordinates using the ~en_distance_rescaled~ and + ~en_distance_rescaled_deriv_e~ which are already calculated previously. + + TODO: write equations. + +*** Get + #+begin_src c :comments org :tangle (eval h_func) :noweb yes +qmckl_exit_code qmckl_get_jastrow_factor_en_deriv_e(qmckl_context context, double* const factor_en_deriv_e); + #+end_src + + #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none +qmckl_exit_code qmckl_get_jastrow_factor_en_deriv_e(qmckl_context context, double* const factor_en_deriv_e) +{ + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { + return QMCKL_NULL_CONTEXT; + } + + qmckl_exit_code rc; + + rc = qmckl_provide_factor_en_deriv_e(context); + if (rc != QMCKL_SUCCESS) return rc; + + qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; + assert (ctx != NULL); + + memcpy(factor_en_deriv_e, ctx->jastrow.factor_en_deriv_e, ctx->electron.walk_num*sizeof(double)); + + return QMCKL_SUCCESS; +} + #+end_src + +*** Provide :noexport: + #+begin_src c :comments org :tangle (eval h_private_func) :noweb yes :exports none +qmckl_exit_code qmckl_provide_factor_en_deriv_e(qmckl_context context); + #+end_src + + #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none +qmckl_exit_code qmckl_provide_factor_en_deriv_e(qmckl_context context) +{ + + qmckl_exit_code rc; + + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { + return QMCKL_NULL_CONTEXT; + } + + qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; + assert (ctx != NULL); + + /* Check if en rescaled distance is provided */ + rc = qmckl_provide_en_distance_rescaled(context); + if(rc != QMCKL_SUCCESS) return rc; + + /* Compute if necessary */ + if (ctx->date > ctx->jastrow.factor_en_deriv_e_date) { + + /* Allocate array */ + if (ctx->jastrow.factor_en_deriv_e == NULL) { + + qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; + mem_info.size = ctx->electron.walk_num * sizeof(double); + double* factor_en_deriv_e = (double*) qmckl_malloc(context, mem_info); + + if (factor_en_deriv_e == NULL) { + return qmckl_failwith( context, + QMCKL_ALLOCATION_FAILED, + "qmckl_provide_factor_en_deriv_e", + NULL); + } + ctx->jastrow.factor_en_deriv_e = factor_en_deriv_e; + } + + qmckl_exit_code rc = + qmckl_compute_factor_en_deriv_e(context, + ctx->electron.walk_num, + ctx->electron.num, + ctx->nucleus.num, + ctx->jastrow.type_nucl_num, + ctx->jastrow.type_nucl_vector, + ctx->jastrow.aord_num, + ctx->jastrow.aord_vector, + ctx->electron.en_distance_rescaled, + ctx->electron.en_distance_rescaled_deriv_e, + ctx->jastrow.factor_en_deriv_e); + if (rc != QMCKL_SUCCESS) { + return rc; + } + + ctx->jastrow.factor_en_deriv_e_date = ctx->date; + } + + return QMCKL_SUCCESS; +} + #+end_src + +*** Compute + :PROPERTIES: + :Name: qmckl_compute_factor_en_deriv_e + :CRetType: qmckl_exit_code + :FRetType: qmckl_exit_code + :END: + + #+NAME: qmckl_factor_en_deriv_e_args + | qmckl_context | context | in | Global state | + | int64_t | walk_num | in | Number of walkers | + | int64_t | elec_num | in | Number of electrons | + | int64_t | nucl_num | in | Number of nucleii | + | int64_t | type_nucl_num | in | Number of unique nuclei | + | int64_t | type_nucl_vector[type_nucl_num] | in | IDs of unique nucleii | + | int64_t | aord_num | in | Number of coefficients | + | double | aord_vector[aord_num + 1][type_nucl_num] | in | List of coefficients | + | double | en_distance_rescaled[walk_num][nucl_num][elec_num] | in | Electron-nucleus distances | + | double | en_distance_rescaled_deriv_e[walk_num][4][nucl_num][elec_num] | in | Electron-nucleus distance derivatives | + | double | factor_en_deriv_e[walk_num][4][elec_num] | out | Electron-nucleus jastrow | + + #+begin_src f90 :comments org :tangle (eval f) :noweb yes +integer function qmckl_compute_factor_en_deriv_e_f(context, walk_num, elec_num, nucl_num, type_nucl_num, & + type_nucl_vector, aord_num, aord_vector, & + en_distance_rescaled, en_distance_rescaled_deriv_e, factor_en_deriv_e) & + result(info) + use qmckl + implicit none + integer(qmckl_context), intent(in) :: context + integer*8 , intent(in) :: walk_num, elec_num, aord_num, nucl_num, type_nucl_num + integer*8 , intent(in) :: type_nucl_vector(type_nucl_num) + double precision , intent(in) :: aord_vector(aord_num, nucl_num) + double precision , intent(in) :: en_distance_rescaled(walk_num, elec_num, nucl_num) + double precision , intent(in) :: en_distance_rescaled_deriv_e(walk_num, 4, elec_num, nucl_num) + double precision , intent(out) :: factor_en_deriv_e(walk_num,4,elec_num) + + integer*8 :: i, a, p, ipar, nw, ii + double precision :: x, spin_fact, den, invden, invden2, invden3, xinv + double precision :: y, lap1, lap2, lap3, third + double precision, dimension(3) :: power_ser_g + double precision, dimension(4) :: dx + + info = QMCKL_SUCCESS + + if (context == QMCKL_NULL_CONTEXT) then + info = QMCKL_INVALID_CONTEXT + return + endif + + if (walk_num <= 0) then + info = QMCKL_INVALID_ARG_2 + return + endif + + if (elec_num <= 0) then + info = QMCKL_INVALID_ARG_3 + return + endif + + if (nucl_num <= 0) then + info = QMCKL_INVALID_ARG_4 + return + endif + + if (aord_num <= 0) then + info = QMCKL_INVALID_ARG_7 + return + endif + + factor_en_deriv_e = 0.0d0 + third = 1.0d0 / 3.0d0 + + do nw =1, walk_num + do a = 1, nucl_num + do i = 1, elec_num + x = en_distance_rescaled(nw, i, a) + if(abs(x) < 1.0d-18) continue + power_ser_g = 0.0d0 + den = 1.0d0 + aord_vector(2, type_nucl_vector(a)) * x + invden = 1.0d0 / den + invden2 = invden * invden + invden3 = invden2 * invden + xinv = 1.0d0 / x + + do ii = 1, 4 + dx(ii) = en_distance_rescaled_deriv_e(nw, ii, i, a) + end do + + lap1 = 0.0d0 + lap2 = 0.0d0 + lap3 = 0.0d0 + do ii = 1, 3 + x = en_distance_rescaled(nw, i, a) + do p = 2, aord_num + y = p * aord_vector(p + 1, type_nucl_vector(a)) * x + power_ser_g(ii) = power_ser_g(ii) + y * dx(ii) + lap1 = lap1 + (p - 1) * y * xinv * dx(ii) * dx(ii) + lap2 = lap2 + y + x = x * en_distance_rescaled(nw, i, a) + end do + + lap3 = lap3 - 2.0d0 * aord_vector(2, type_nucl_vector(a)) * dx(ii) * dx(ii) + + factor_en_deriv_e(nw, ii, i) = factor_en_deriv_e(nw, ii, i) + aord_vector(1, type_nucl_vector(a)) & + * dx(ii) * invden2 & + + power_ser_g(ii) + + end do + + ii = 4 + lap2 = lap2 * dx(ii) * third + lap3 = lap3 + den * dx(ii) + lap3 = lap3 * aord_vector(1, type_nucl_vector(a)) * invden3 + factor_en_deriv_e(nw, ii, i) = factor_en_deriv_e(nw, ii, i) + lap1 + lap2 + lap3 + + end do + end do + end do + +end function qmckl_compute_factor_en_deriv_e_f + #+end_src + + #+CALL: generate_c_header(table=qmckl_factor_en_deriv_e_args,rettyp=get_value("CRetType"),fname=get_value("Name")) + + #+RESULTS: + #+begin_src c :tangle (eval h_func) :comments org + qmckl_exit_code qmckl_compute_factor_en_deriv_e ( + 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* aord_vector, + const double* en_distance_rescaled, + const double* en_distance_rescaled_deriv_e, + double* const factor_en_deriv_e ); + #+end_src + + + #+CALL: generate_c_interface(table=qmckl_factor_en_deriv_e_args,rettyp=get_value("CRetType"),fname=get_value("Name")) + + #+RESULTS: + #+begin_src f90 :tangle (eval f) :comments org :exports none + integer(c_int32_t) function qmckl_compute_factor_en_deriv_e & + (context, & + walk_num, & + elec_num, & + nucl_num, & + type_nucl_num, & + type_nucl_vector, & + aord_num, & + aord_vector, & + en_distance_rescaled, & + en_distance_rescaled_deriv_e, & + factor_en_deriv_e) & + bind(C) result(info) + + use, intrinsic :: iso_c_binding + implicit none + + integer (c_int64_t) , intent(in) , value :: context + integer (c_int64_t) , intent(in) , value :: walk_num + integer (c_int64_t) , intent(in) , value :: elec_num + integer (c_int64_t) , intent(in) , value :: nucl_num + integer (c_int64_t) , intent(in) , value :: type_nucl_num + integer (c_int64_t) , intent(in) :: type_nucl_vector(type_nucl_num) + integer (c_int64_t) , intent(in) , value :: aord_num + real (c_double ) , intent(in) :: aord_vector(type_nucl_num,aord_num + 1) + real (c_double ) , intent(in) :: en_distance_rescaled(elec_num,nucl_num,walk_num) + real (c_double ) , intent(in) :: en_distance_rescaled_deriv_e(elec_num,nucl_num,4,walk_num) + real (c_double ) , intent(out) :: factor_en_deriv_e(elec_num,4,walk_num) + + integer(c_int32_t), external :: qmckl_compute_factor_en_deriv_e_f + info = qmckl_compute_factor_en_deriv_e_f & + (context, & + walk_num, & + elec_num, & + nucl_num, & + type_nucl_num, & + type_nucl_vector, & + aord_num, & + aord_vector, & + en_distance_rescaled, & + en_distance_rescaled_deriv_e, & + factor_en_deriv_e) + + end function qmckl_compute_factor_en_deriv_e + #+end_src + +*** Test + #+begin_src python :results output :exports none :noweb yes +import numpy as np + +<> + +factor_en = 0.0 +for a in range(0,nucl_num): + for i in range(0,elec_num): + x = en_distance_rescaled[i][a] + pow_ser = 0.0 + + for p in range(2,aord_num+1): + x = x * en_distance_rescaled[i][a] + pow_ser = pow_ser + aord_vector[(p-1) + 1][type_nucl_vector[a]-1] * x + + factor_en = factor_en + aord_vector[0][type_nucl_vector[a]-1] * en_distance_rescaled[i][a] \ + / (1.0 + aord_vector[1][type_nucl_vector[a]-1] * en_distance_rescaled[i][a]) \ + + pow_ser +print("factor_en :",factor_en) + + #+end_src + + #+RESULTS: + : factor_en : -5.865822569188727 + + + #+begin_src c :tangle (eval c_test) +/* Check if Jastrow is properly initialized */ +assert(qmckl_jastrow_provided(context)); + +//double factor_en[walk_num]; +//rc = qmckl_get_jastrow_factor_en(context, factor_en); +// +//// calculate factor_en +//assert(fabs(factor_en[0]+5.865822569188727) < 1.e-12); + + #+end_src + * End of files :noexport: