1
0
mirror of https://github.com/TREX-CoE/qmckl.git synced 2024-06-30 00:44:52 +02:00

Working on factor_een. #22

This commit is contained in:
vijay gopal chilkuri 2021-07-07 11:23:04 +05:30
parent 2af6e36252
commit a7c1fe526e
2 changed files with 361 additions and 49 deletions

View File

@ -63,33 +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][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 |
| ~een_distance_rescaled~ | ~double[walk_num][nucl_num][num]~ | Electron-nucleus distances |
| ~een_distance_rescaled_date~ | ~uint64_t~ | Last modification date of the electron-electron distances |
| ~een_distance_rescaled_deriv_e~ | ~double[walk_num][4][nucl_num][num]~ | Electron-electron rescaled distances derivatives |
| ~een_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
@ -104,12 +100,10 @@ typedef struct qmckl_electron_struct {
int64_t coord_new_date;
int64_t ee_distance_date;
int64_t en_distance_date;
int64_t een_distance_date;
int64_t ee_distance_rescaled_date;
int64_t ee_distance_rescaled_deriv_e_date;
int64_t en_distance_rescaled_date;
int64_t en_distance_rescaled_deriv_e_date;
int64_t een_distance_rescaled_deriv_e_date;
double* coord_new;
double* coord_old;
double* ee_distance;
@ -118,8 +112,6 @@ typedef struct qmckl_electron_struct {
double* ee_distance_rescaled_deriv_e;
double* en_distance_rescaled;
double* en_distance_rescaled_deriv_e;
double* een_distance_rescaled;
double* een_distance_rescaled_deriv_e;
int32_t uninitialized;
bool provided;
} qmckl_electron_struct;
@ -1559,7 +1551,6 @@ rc = qmckl_get_electron_ee_distance_rescaled_deriv_e(context, ee_distance_rescal
#+end_src
** Electron-nucleus distances
*** Get

View File

@ -83,28 +83,32 @@ int main() {
| ~double~ | ~bord_vector[bord_num + 1]~ | in | Order of b polynomial coefficients |
| ~double~ | ~cord_vector[cord_num][type_nucl_num]~ | in | Order of c polynomial coefficients |
| ~double~ | ~factor_ee[walk_num]~ | out | Jastrow factor: electron-electron part |
| ~double~ | ~factor_ee_date~ | out | Jastrow factor: electron-electron part |
| ~int64_t~ | ~factor_ee_date~ | out | Jastrow factor: electron-electron part |
| ~double~ | ~factor_en[walk_num]~ | out | Jastrow factor: electron-nucleus part |
| ~double~ | ~factor_en_date~ | out | Jastrow factor: electron-nucleus part |
| ~int64_t~ | ~factor_en_date~ | out | Jastrow factor: electron-nucleus part |
| ~double~ | ~factor_een[walk_num]~ | out | Jastrow factor: electron-electron-nucleus part |
| ~double~ | ~factor_een_date~ | out | Jastrow factor: electron-electron-nucleus part |
| ~int64_t~ | ~factor_een_date~ | out | Jastrow factor: electron-electron-nucleus part |
| ~double~ | ~factor_ee_deriv_e[4][nelec][walk_num]~ | out | Derivative of the Jastrow factor: electron-electron-nucleus part |
| ~double~ | ~factor_ee_deriv_e_date~ | out | Keep track of the date for the derivative |
| ~int64_t~ | ~factor_ee_deriv_e_date~ | out | Keep track of the date for the derivative |
| ~double~ | ~factor_en_deriv_e[4][nelec][walk_num]~ | out | Derivative of the Jastrow factor: electron-electron-nucleus part |
| ~double~ | ~factor_en_deriv_e_date~ | out | Keep track of the date for the en derivative |
| ~int64_t~ | ~factor_en_deriv_e_date~ | out | Keep track of the date for the en derivative |
| ~double~ | ~factor_een_deriv_e[4][nelec][walk_num]~ | out | Derivative of the Jastrow factor: electron-electron-nucleus part |
| ~double~ | ~factor_een_deriv_e_date~ | out | Keep track of the date for the een derivative |
| ~int64_t~ | ~factor_een_deriv_e_date~ | out | Keep track of the date for the een derivative |
computed data:
|-----------+-------------------------------------------------------------+-------------------------------------------------|
| ~int64_t~ | ~dim_cord_vec~ | Number of unique C coefficients |
| ~double~ | ~asymp_jasb[2]~ | Asymptotic component |
| ~int64_t~ | ~asymp_jasb_date~ | Asymptotic component |
| ~double~ | ~coord_vect_full[dim_cord_vec][nucl_num]~ | vector of non-zero coefficients |
| ~int64_t~ | ~lkpm_of_cindex[4][dim_cord_vec]~ | Transform l,k,p, and m into consecutive indices |
| ~double~ | ~tmp_c[elec_num][nucl_num][ncord + 1][ncord][walk_num]~ | vector of non-zero coefficients |
| ~double~ | ~dtmp_c[elec_num][4][nucl_num][ncord + 1][ncord][walk_num]~ | vector of non-zero coefficients |
|-----------+-------------------------------------------------------------+-------------------------------------------------------------------------------|
| ~int64_t~ | ~dim_cord_vec~ | Number of unique C coefficients |
| ~double~ | ~asymp_jasb[2]~ | Asymptotic component |
| ~int64_t~ | ~asymp_jasb_date~ | Asymptotic component |
| ~double~ | ~coord_vect_full[dim_cord_vec][nucl_num]~ | vector of non-zero coefficients |
| ~int64_t~ | ~lkpm_of_cindex[4][dim_cord_vec]~ | Transform l,k,p, and m into consecutive indices |
| ~double~ | ~tmp_c[elec_num][nucl_num][ncord + 1][ncord][walk_num]~ | vector of non-zero coefficients |
| ~double~ | ~dtmp_c[elec_num][4][nucl_num][ncord + 1][ncord][walk_num]~ | vector of non-zero coefficients |
| ~double~ | ~een_rescaled_e[walk_num][elec_num][elec_num][0:cord_num]~ | The electron-electron rescaled distances raised to the powers defined by cord |
| ~int64_t~ | ~een_rescaled_e_date~ | Keep track of the date of creation |
| ~double~ | ~een_rescaled_n[walk_num][elec_num][nucl_num][0:cord_num]~ | The electron-electron rescaled distances raised to the powers defined by cord |
| ~int64_t~ | ~een_rescaled_n_date~ | Keep track of the date of creation |
For H2O we have the following data:
@ -293,6 +297,10 @@ typedef struct qmckl_jastrow_struct{
double * coord_vect_full;
double * tmp_c;
double * dtmp_c;
double * een_rescaled_e;
double * een_rescaled_n;
int64_t een_rescaled_e_date;
int64_t een_rescaled_n_date;
bool provided;
char * type;
} qmckl_jastrow_struct;
@ -1147,7 +1155,7 @@ assert(qmckl_nucleus_provided(context));
the dependencies are more recent than the date of the data to
compute. If it is the case, then the data is recomputed and the
current date is stored.
** Asymptotic component for \(f_{ee}\)
Calculate the asymptotic component ~asymp_jasb~ to be substracted from the final
@ -1411,9 +1419,8 @@ assert(fabs(asymp_jasb[1]-0.31567342786262853) < 1.e-12);
#+end_src
** Electron-electron component \(f_{ee}\)
Calculate the electron-electron jastrow component ~factor_ee~ using the ~asymp_jasb~
componenet and the electron-electron rescaled distances ~ee_distance_rescaled~.
@ -2829,6 +2836,320 @@ assert(fabs(factor_en_deriv_e[0][3][0]+0.9667363412285741 ) < 1.e-12);
#+end_src
** Electron-electron rescaled distances for each order
~een_rescaled_e~ stores the table of the rescaled distances between all
pairs of electrons and raised to the power \(p\) defined by ~cord_num~:
\[
C_{ij,p} = \left( 1 - \exp{-\kappa C_{ij}} \right)^p
\]
where \(C_{ij}\) is the matrix of electron-electron distances.
*** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_een_rescaled_e(qmckl_context context, double* const distance_rescaled);
#+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_een_rescaled_e(qmckl_context context, double* const distance_rescaled)
{
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT;
}
qmckl_exit_code rc;
rc = qmckl_provide_een_rescaled_e(context);
if (rc != QMCKL_SUCCESS) return rc;
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
size_t sze = ctx->electron.num * ctx->electron.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1);
memcpy(distance_rescaled, ctx->jastrow.een_rescaled_e, sze * 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_een_rescaled_e(qmckl_context context);
#+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_provide_een_rescaled_e(qmckl_context context)
{
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 ee distance is provided */
qmckl_exit_code rc = qmckl_provide_ee_distance(context);
if(rc != QMCKL_SUCCESS) return rc;
/* Compute if necessary */
if (ctx->date > ctx->jastrow.een_rescaled_e_date) {
/* Allocate array */
if (ctx->jastrow.een_rescaled_e == NULL) {
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
mem_info.size = ctx->electron.num * ctx->electron.num *
ctx->electron.walk_num * (ctx->jastrow.cord_num + 1) * sizeof(double);
double* een_rescaled_e = (double*) qmckl_malloc(context, mem_info);
if (een_rescaled_e == NULL) {
return qmckl_failwith( context,
QMCKL_ALLOCATION_FAILED,
"qmckl_een_rescaled_e",
NULL);
}
ctx->jastrow.een_rescaled_e = een_rescaled_e;
}
qmckl_exit_code rc =
qmckl_compute_een_rescaled_e(context,
ctx->electron.walk_num,
ctx->electron.num,
ctx->jastrow.cord_num,
ctx->electron.rescale_factor_kappa_ee,
ctx->electron.ee_distance,
ctx->jastrow.een_rescaled_e);
if (rc != QMCKL_SUCCESS) {
return rc;
}
ctx->jastrow.een_rescaled_e_date = ctx->date;
}
return QMCKL_SUCCESS;
}
#+end_src
*** Compute
:PROPERTIES:
:Name: qmckl_compute_een_rescaled_e
:CRetType: qmckl_exit_code
:FRetType: qmckl_exit_code
:END:
#+NAME: qmckl_factor_een_rescaled_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 | cord_num | in | Order of polynomials |
| double | rescale_factor_kappa_ee | in | Factor to rescale ee distances |
| double | ee_distance[walk_num][elec_num][elec_num] | in | Electron-electron distances |
| double | een_rescaled_e[walk_num][elec_num][elec_num][0:cord_num] | out | Electron-electron rescaled distances |
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
integer function qmckl_compute_een_rescaled_e_f(context, walk_num, elec_num, cord_num, rescale_factor_kappa_ee, &
ee_distance, een_rescaled_e) &
result(info)
use qmckl
implicit none
integer(qmckl_context), intent(in) :: context
integer*8 , intent(in) :: walk_num
integer*8 , intent(in) :: elec_num
integer*8 , intent(in) :: cord_num
double precision , intent(in) :: rescale_factor_kappa_ee
double precision , intent(in) :: ee_distance(elec_num,elec_num,walk_num)
double precision , intent(out) :: een_rescaled_e(0:cord_num,elec_num,elec_num,walk_num)
double precision,dimension(:,:),allocatable :: een_rescaled_e_ij
double precision :: x
integer*8 :: i, j, k, l, nw
allocate(een_rescaled_e_ij(elec_num * (elec_num - 1) / 2, cord_num + 1))
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 (cord_num <= 0) then
info = QMCKL_INVALID_ARG_4
return
endif
een_rescaled_e_ij(:, 1) = 1.0d0
! Prepare table of exponentiated distances raised to appropriate power
do nw = 1, walk_num
k = 0
do j = 1, elec_num
do i = 1, j - 1
k = k + 1
een_rescaled_e_ij(k, 2) = dexp(-rescale_factor_kappa_ee * ee_distance(i, j, nw))
end do
end do
do l = 2, cord_num
do k = 1, elec_num * (elec_num - 1)/2
een_rescaled_e_ij(k, l + 1) = een_rescaled_e_ij(k, l + 1 - 1) * een_rescaled_e_ij(k, 1)
end do
end do
! prepare the actual een table
een_rescaled_e = 0.0d0
een_rescaled_e(0, :, :, :) = 1.0d0
do l = 1, cord_num
k = 0
do j = 1, elec_num
do i = 1, j - 1
k = k + 1
x = een_rescaled_e_ij(k, l + 1)
een_rescaled_e(l, i, j, nw) = x
een_rescaled_e(l, j, i, nw) = x
end do
end do
end do
end do
end function qmckl_compute_een_rescaled_e_f
#+end_src
#+CALL: generate_c_header(table=qmckl_factor_een_rescaled_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_een_rescaled_e (
const qmckl_context context,
const int64_t walk_num,
const int64_t elec_num,
const int64_t cord_num,
const double rescale_factor_kappa_ee,
const double* ee_distance,
double* const een_rescaled_e );
#+end_src
#+CALL: generate_c_interface(table=qmckl_factor_een_rescaled_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_een_rescaled_e &
(context, walk_num, elec_num, cord_num, rescale_factor_kappa_ee, ee_distance, een_rescaled_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 :: cord_num
real (c_double ) , intent(in) , value :: rescale_factor_kappa_ee
real (c_double ) , intent(in) :: ee_distance(elec_num,elec_num,walk_num)
real (c_double ) , intent(out) :: een_rescaled_e(0:cord_num,elec_num,elec_num,walk_num)
integer(c_int32_t), external :: qmckl_compute_een_rescaled_e_f
info = qmckl_compute_een_rescaled_e_f &
(context, walk_num, elec_num, cord_num, rescale_factor_kappa_ee, ee_distance, een_rescaled_e)
end function qmckl_compute_een_rescaled_e
#+end_src
*** Test
#+begin_src python :results output :exports none :noweb yes
import numpy as np
<<jastrow_data>>
elec_coord = np.array(elec_coord)[0]
elec_dist = np.zeros(shape=(elec_num, elec_num),dtype=float)
for i in range(elec_num):
for j in range(elec_num):
elec_dist[i, j] = np.linalg.norm(elec_coord[i] - elec_coord[j])
kappa = 1.0
een_rescaled_e_ij = np.zeros(shape=(elec_num * (elec_num - 1)//2, cord_num+1), dtype=float)
een_rescaled_e_ij[:,0] = 1.0
k = 0
for j in range(elec_num):
for i in range(j - 1):
een_rescaled_e_ij[k, 1] = np.exp(-kappa * elec_dist[i, j])
k = k + 1
for l in range(2, cord_num + 1):
for k in range(elec_num * (elec_num - 1)//2):
een_rescaled_e_ij[k, l] = een_rescaled_e_ij[k, l - 1] * een_rescaled_e_ij[k, 1]
een_rescaled_e = np.zeros(shape=(elec_num, elec_num, cord_num + 1), dtype=float)
een_rescaled_e[:,:,0] = 1.0
for l in range(1,cord_num+1):
k = 0
for j in range(elec_num):
for i in range(j - 1):
x = een_rescaled_e_ij[k, l]
een_rescaled_e[i, j, l] = x
een_rescaled_e[j, i, l] = x
k = k + 1
print(" een_rescaled_e[0, 2, 1] = ",een_rescaled_e[0, 2, 1])
print(" een_rescaled_e[0, 3, 1] = ",een_rescaled_e[0, 3, 1])
print(" een_rescaled_e[0, 4, 1] = ",een_rescaled_e[0, 4, 1])
print(" een_rescaled_e[1, 3, 2] = ",een_rescaled_e[1, 3, 2])
print(" een_rescaled_e[1, 4, 2] = ",een_rescaled_e[1, 4, 2])
print(" een_rescaled_e[1, 5, 2] = ",een_rescaled_e[1, 5, 2])
#+end_src
#+RESULTS:
: een_rescaled_e[0, 2, 1] = 0.08084493981483197
: een_rescaled_e[0, 3, 1] = 0.1066745707571846
: een_rescaled_e[0, 4, 1] = 0.01754273169464735
: een_rescaled_e[1, 3, 2] = 0.02214680362033448
: een_rescaled_e[1, 4, 2] = 0.0005700154999202759
: een_rescaled_e[1, 5, 2] = 0.3424402276009091
#+begin_src c :tangle (eval c_test)
assert(qmckl_electron_provided(context));
double een_rescaled_e[walk_num][elec_num][elec_num][(cord_num + 1)];
rc = qmckl_get_jastrow_een_rescaled_e(context, &(een_rescaled_e[0][0][0][0]));
// value of (0,2,1)
//printf("%10.15f = \n", een_rescaled_e[0][0][2][1]);
//printf("%10.15f = \n", een_rescaled_e[0][0][3][1]);
//printf("%10.15f = \n", een_rescaled_e[0][0][4][1]);
//printf("%10.15f = \n", een_rescaled_e[0][1][3][2]);
//printf("%10.15f = \n", een_rescaled_e[0][1][4][2]);
//printf("%10.15f = \n", een_rescaled_e[0][1][5][2]);
//assert(fabs(een_rescaled_e[0][0][2][1]-) < 1.e-12);
//assert(fabs(een_rescaled_e[0][0][3][1]-) < 1.e-12);
//assert(fabs(een_rescaled_e[0][0][4][1]-) < 1.e-12);
//assert(fabs(een_rescaled_e[0][1][3][2]-) < 1.e-12);
//assert(fabs(een_rescaled_e[0][1][4][2]-) < 1.e-12);
//assert(fabs(een_rescaled_e[0][1][5][2]-) < 1.e-12);
#+end_src
* End of files :noexport:
#+begin_src c :tangle (eval h_private_type)