TREX
/
qmckl
mirror of https://github.com/TREX-CoE/qmckl.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity

Transposed en_distance

This commit is contained in:
Anthony Scemama 2023-03-13 15:32:35 +01:00
parent 44c4c6c6d5
commit 3f33db6887
4 changed files with 90 additions and 95 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
org/qmckl_ao.org
View File

@ -6305,7 +6305,7 @@ assert( fabs(ao_value[26][224] - ( 7.175045873560788e-10)) < 1.e-14 );
:CRetType: qmckl_exit_code
:FRetType: qmckl_exit_code
:END:
*** Unoptimized version
*** Reference version
#+NAME: qmckl_ao_vgl_args_doc
| Variable | Type | In/Out | Description |
|-----------------------+-----------------------------------+--------+----------------------------------------------|

43
org/qmckl_electron.org
View File

@ -95,7 +95,7 @@ int main() {
|-------------------------------------+----------------------------------------+----------------------------------------------------------------------|
| ~ee_distance~ | ~double[walker.num][num][num]~ | Electron-electron distances |
| ~ee_distance_date~ | ~uint64_t~ | Last modification date of the electron-electron distances |
| ~en_distance~ | ~double[walker.num][nucl_num][num]~ | Electron-nucleus distances |
| ~en_distance~ | ~double[walker.num][num][nucl_num]~ | Electron-nucleus distances |
| ~en_distance_date~ | ~uint64_t~ | Last modification date of the electron-electron distances |
| ~ee_potential~ | ~double[walker.num]~ | Electron-electron potential energy |
| ~ee_potential_date~ | ~uint64_t~ | Last modification date of the electron-electron potential |
@ -1235,7 +1235,7 @@ qmckl_exit_code qmckl_provide_en_distance(qmckl_context context)
| ~walk_num~ | ~int64_t~ | in | Number of walkers |
| ~elec_coord~ | ~double[3][walk_num][elec_num]~ | in | Electron coordinates |
| ~nucl_coord~ | ~double[3][elec_num]~ | in | Nuclear coordinates |
| ~en_distance~ | ~double[walk_num][nucl_num][elec_num]~ | out | Electron-nucleus distances |
| ~en_distance~ | ~double[walk_num][elec_num][nucl_num]~ | out | Electron-nucleus distances |
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
integer function qmckl_compute_en_distance_f(context, elec_num, nucl_num, walk_num, elec_coord, nucl_coord, en_distance) &
@ -1248,7 +1248,7 @@ integer function qmckl_compute_en_distance_f(context, elec_num, nucl_num, walk_n
integer*8 , intent(in) :: walk_num
double precision , intent(in) :: elec_coord(elec_num,walk_num,3)
double precision , intent(in) :: nucl_coord(nucl_num,3)
double precision , intent(out) :: en_distance(elec_num,nucl_num,walk_num)
double precision , intent(out) :: en_distance(nucl_num,elec_num,walk_num)
integer*8 :: k
@ -1274,15 +1274,10 @@ integer function qmckl_compute_en_distance_f(context, elec_num, nucl_num, walk_n
return
endif
do k=1,walk_num
info = qmckl_distance(context, 'T', 'T', elec_num, nucl_num, &
elec_coord(1,k,1), elec_num * walk_num, &
info = qmckl_distance(context, 'T', 'T', nucl_num, elec_num * walk_num, &
nucl_coord, nucl_num, &
en_distance(1,1,k), elec_num)
if (info /= QMCKL_SUCCESS) then
exit
endif
end do
elec_coord, elec_num * walk_num, &
en_distance, nucl_num)
end function qmckl_compute_en_distance_f
#+end_src
@ -1315,7 +1310,7 @@ qmckl_exit_code qmckl_compute_en_distance (
integer (c_int64_t) , intent(in) , value :: walk_num
real (c_double ) , intent(in) :: elec_coord(elec_num,walk_num,3)
real (c_double ) , intent(in) :: nucl_coord(elec_num,3)
real (c_double ) , intent(out) :: en_distance(elec_num,nucl_num,walk_num)
real (c_double ) , intent(out) :: en_distance(nucl_num,elec_num,walk_num)
integer(c_int32_t), external :: qmckl_compute_en_distance_f
info = qmckl_compute_en_distance_f &
@ -1368,7 +1363,7 @@ qmckl_check(context, rc);
assert(qmckl_nucleus_provided(context));
double en_distance[walk_num][nucl_num][elec_num];
double en_distance[walk_num][elec_num][nucl_num];
rc = qmckl_get_electron_en_distance(context, &(en_distance[0][0][0]));
qmckl_check(context, rc);
@ -1378,19 +1373,19 @@ qmckl_check(context, rc);
assert(fabs(en_distance[0][0][0] - 7.546738741619978) < 1.e-12);
// (1,2,1)
assert(fabs(en_distance[0][1][0] - 8.77102435246984) < 1.e-12);
assert(fabs(en_distance[0][0][1] - 8.77102435246984) < 1.e-12);
// (2,1,1)
assert(fabs(en_distance[0][0][1] - 3.698922010513608) < 1.e-12);
assert(fabs(en_distance[0][1][0] - 3.698922010513608) < 1.e-12);
// (1,1,2)
assert(fabs(en_distance[1][0][0] - 5.824059436060509) < 1.e-12);
// (1,2,2)
assert(fabs(en_distance[1][1][0] - 7.080482110317645) < 1.e-12);
assert(fabs(en_distance[1][0][1] - 7.080482110317645) < 1.e-12);
// (2,1,2)
assert(fabs(en_distance[1][0][1] - 3.1804527583077356) < 1.e-12);
assert(fabs(en_distance[1][1][0] - 3.1804527583077356) < 1.e-12);
#+end_src
@ -1512,7 +1507,7 @@ qmckl_exit_code qmckl_provide_en_potential(qmckl_context context)
| ~nucl_num~ | ~int64_t~ | in | Number of nuclei |
| ~walk_num~ | ~int64_t~ | in | Number of walkers |
| ~charge~ | ~double[nucl_num]~ | in | charge of nucleus |
| ~en_distance~ | ~double[walk_num][nucl_num][elec_num]~ | in | Electron-electron rescaled distances |
| ~en_distance~ | ~double[walk_num][elec_num][nucl_num]~ | in | Electron-electron distances |
| ~en_potential~ | ~double[walk_num]~ | out | Electron-electron potential |
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
@ -1526,7 +1521,7 @@ integer function qmckl_compute_en_potential_f(context, elec_num, nucl_num, walk_
integer*8 , intent(in) :: nucl_num
integer*8 , intent(in) :: walk_num
double precision , intent(in) :: charge(nucl_num)
double precision , intent(in) :: en_distance(elec_num,nucl_num,walk_num)
double precision , intent(in) :: en_distance(nucl_num,elec_num,walk_num)
double precision , intent(out) :: en_potential(walk_num)
integer*8 :: nw, i, j
@ -1550,10 +1545,10 @@ integer function qmckl_compute_en_potential_f(context, elec_num, nucl_num, walk_
en_potential = 0.0d0
do nw=1,walk_num
do j=1,nucl_num
do i=1,elec_num
if (dabs(en_distance(i,j,nw)) > 1e-5) then
en_potential(nw) = en_potential(nw) - charge(j)/(en_distance(i,j,nw))
do i=1,elec_num
do j=1,nucl_num
if (dabs(en_distance(j,i,nw)) > 1.d-6) then
en_potential(nw) = en_potential(nw) - charge(j)/(en_distance(j,i,nw))
endif
end do
end do
@ -1592,7 +1587,7 @@ end function qmckl_compute_en_potential_f
integer (c_int64_t) , intent(in) , value :: nucl_num
integer (c_int64_t) , intent(in) , value :: walk_num
real (c_double ) , intent(in) :: charge(nucl_num)
real (c_double ) , intent(in) :: en_distance(elec_num,nucl_num,walk_num)
real (c_double ) , intent(in) :: en_distance(nucl_num,elec_num,walk_num)
real (c_double ) , intent(out) :: en_potential(walk_num)
integer(c_int32_t), external :: qmckl_compute_en_potential_f

134
org/qmckl_jastrow.org
View File

@ -11,7 +11,7 @@
\[
J(\mathbf{r},\mathbf{R}) = J_{\text{eN}}(\mathbf{r},\mathbf{R}) + J_{\text{ee}}(\mathbf{r}) + J_{\text{eeN}}(\mathbf{r},\mathbf{R})
\]
In the following, we use the notations $r_{ij} = |\mathbf{r}_i - \mathbf{r}_j|$ and
$R_{i\alpha} = |\mathbf{r}_i - \mathbf{R}_\alpha|$.
@ -19,7 +19,7 @@
\[
J_{\text{eN}}(\mathbf{r},\mathbf{R}) =
\sum_{\alpha=1}^{N_\text{nucl}} \sum_{i=1}^{N_\text{elec}}
\sum_{\alpha=1}^{N_\text{nucl}} \sum_{i=1}^{N_\text{elec}}
\frac{a_{1\,\alpha}\, f_\alpha(R_{i\,\alpha})}{1+a_{2\,\alpha}\, f_\alpha(R_{i\alpha})} +
\sum_{p=2}^{N_\text{ord}^a} a_{p+1\,\alpha}\, [f_\alpha(R_{i\alpha})]^p - J_{eN}^{\infty \alpha}
\]
@ -640,7 +640,7 @@ qmckl_set_jastrow_a_vector(qmckl_context context,
"qmckl_set_jastrow_a_vector",
"aord_num not initialized");
}
int64_t type_nucl_num = ctx->jastrow.type_nucl_num;
if (type_nucl_num <= 0) {
@ -825,7 +825,7 @@ qmckl_set_jastrow_rescale_factor_ee(qmckl_context context,
const double rescale_factor_ee) {
int32_t mask = 1 << 8;
<<pre2>>
if (rescale_factor_ee <= 0.0) {
@ -979,7 +979,7 @@ interface
integer(c_int64_t), intent(in), value :: size_max
double precision, intent(in) :: kappa_en(size_max)
end function qmckl_set_jastrow_rescale_factor_en
integer(qmckl_exit_code) function qmckl_set_jastrow_aord_num (context, &
aord_num) bind(C)
use, intrinsic :: iso_c_binding
@ -988,7 +988,7 @@ interface
integer (qmckl_context) , intent(in) , value :: context
integer(c_int64_t), intent(in), value :: aord_num
end function qmckl_set_jastrow_aord_num
integer(qmckl_exit_code) function qmckl_set_jastrow_bord_num (context, &
bord_num) bind(C)
use, intrinsic :: iso_c_binding
@ -997,7 +997,7 @@ interface
integer (qmckl_context) , intent(in) , value :: context
integer(c_int64_t), intent(in), value :: bord_num
end function qmckl_set_jastrow_bord_num
integer(qmckl_exit_code) function qmckl_set_jastrow_cord_num (context, &
cord_num) bind(C)
use, intrinsic :: iso_c_binding
@ -1015,7 +1015,7 @@ interface
integer (qmckl_context) , intent(in) , value :: context
integer(c_int64_t), intent(in), value :: type_nucl_num
end function qmckl_set_jastrow_type_nucl_num
integer(qmckl_exit_code) function qmckl_set_jastrow_type_nucl_vector (context, &
type_nucl_vector, size_max) bind(C)
use, intrinsic :: iso_c_binding
@ -1025,7 +1025,7 @@ interface
integer(c_int64_t), intent(in), value :: size_max
integer(c_int64_t), intent(in) :: type_nucl_vector(size_max)
end function qmckl_set_jastrow_type_nucl_vector
integer(qmckl_exit_code) function qmckl_set_jastrow_a_vector(context, &
a_vector, size_max) bind(C)
use, intrinsic :: iso_c_binding
@ -1035,7 +1035,7 @@ interface
integer(c_int64_t), intent(in), value :: size_max
double precision, intent(in) :: a_vector(size_max)
end function qmckl_set_jastrow_a_vector
integer(qmckl_exit_code) function qmckl_set_jastrow_b_vector(context, &
b_vector, size_max) bind(C)
use, intrinsic :: iso_c_binding
@ -1045,7 +1045,7 @@ interface
integer(c_int64_t), intent(in), value :: size_max
double precision, intent(in) :: b_vector(size_max)
end function qmckl_set_jastrow_b_vector
integer(qmckl_exit_code) function qmckl_set_jastrow_c_vector(context, &
c_vector, size_max) bind(C)
use, intrinsic :: iso_c_binding
@ -1055,7 +1055,7 @@ interface
integer(c_int64_t), intent(in), value :: size_max
double precision, intent(in) :: c_vector(size_max)
end function qmckl_set_jastrow_c_vector
end interface
#+end_src
@ -1074,7 +1074,7 @@ qmckl_exit_code qmckl_get_jastrow_rescale_factor_ee (const qmckl_context contex
qmckl_exit_code qmckl_get_jastrow_rescale_factor_en (const qmckl_context context, double* const rescale_factor_en, const int64_t size_max);
#+end_src
Along with these core functions, calculation of the jastrow factor
requires the following additional information to be set:
@ -1455,7 +1455,7 @@ interface
integer(c_int64_t), intent(in), value :: size_max
double precision, intent(out) :: kappa_en(size_max)
end function qmckl_get_jastrow_rescale_factor_en
integer(qmckl_exit_code) function qmckl_get_jastrow_aord_num (context, &
aord_num) bind(C)
use, intrinsic :: iso_c_binding
@ -1464,7 +1464,7 @@ interface
integer (qmckl_context) , intent(in), value :: context
integer(c_int64_t), intent(out) :: aord_num
end function qmckl_get_jastrow_aord_num
integer(qmckl_exit_code) function qmckl_get_jastrow_bord_num (context, &
bord_num) bind(C)
use, intrinsic :: iso_c_binding
@ -1473,7 +1473,7 @@ interface
integer (qmckl_context) , intent(in), value :: context
integer(c_int64_t), intent(out) :: bord_num
end function qmckl_get_jastrow_bord_num
integer(qmckl_exit_code) function qmckl_get_jastrow_cord_num (context, &
cord_num) bind(C)
use, intrinsic :: iso_c_binding
@ -1491,7 +1491,7 @@ interface
integer (qmckl_context) , intent(in), value :: context
integer(c_int64_t), intent(out) :: type_nucl_num
end function qmckl_get_jastrow_type_nucl_num
integer(qmckl_exit_code) function qmckl_get_jastrow_type_nucl_vector (context, &
type_nucl_vector, size_max) bind(C)
use, intrinsic :: iso_c_binding
@ -1501,7 +1501,7 @@ interface
integer(c_int64_t), intent(in), value :: size_max
integer(c_int64_t), intent(out) :: type_nucl_vector(size_max)
end function qmckl_get_jastrow_type_nucl_vector
integer(qmckl_exit_code) function qmckl_get_jastrow_a_vector(context, &
a_vector, size_max) bind(C)
use, intrinsic :: iso_c_binding
@ -1511,7 +1511,7 @@ interface
integer(c_int64_t), intent(in), value :: size_max
double precision, intent(out) :: a_vector(size_max)
end function qmckl_get_jastrow_a_vector
integer(qmckl_exit_code) function qmckl_get_jastrow_b_vector(context, &
b_vector, size_max) bind(C)
use, intrinsic :: iso_c_binding
@ -1521,7 +1521,7 @@ interface
integer(c_int64_t), intent(in), value :: size_max
double precision, intent(out) :: b_vector(size_max)
end function qmckl_get_jastrow_b_vector
integer(qmckl_exit_code) function qmckl_get_jastrow_c_vector(context, &
c_vector, size_max) bind(C)
use, intrinsic :: iso_c_binding
@ -1531,7 +1531,7 @@ interface
integer(c_int64_t), intent(in), value :: size_max
double precision, intent(out) :: c_vector(size_max)
end function qmckl_get_jastrow_c_vector
end interface
#+end_src
@ -1874,7 +1874,7 @@ qmckl_exit_code qmckl_compute_jastrow_asymp_jasb_doc (const qmckl_context contex
const double rescale_factor_ee,
double* const asymp_jasb );
#+end_src
#+begin_src c :comments org :tangle (eval h_private_func) :noweb yes :exports none
qmckl_exit_code qmckl_compute_jastrow_asymp_jasb_hpc (const qmckl_context context,
@ -2067,7 +2067,7 @@ assert(fabs(asymp_jasb[1]-0.31567342786262853) < 1.e-12);
\[
f_\text{ee} = \sum_{i,j<i} \left[
\frac{\delta_{ij}^{\uparrow\downarrow} B_0\, C_{ij}}{1 - B_1\,
C_{ij}} + \sum_{k=2}^{n_\text{ord}} B_k\, C_{ij}^k - {J_{\text{ee}}^{\infty}}_{ij} \right]
C_{ij}} + \sum_{k=2}^{n_\text{ord}} B_k\, C_{ij}^k - {J_{\text{ee}}^{\infty}}_{ij} \right]
\]
$\delta$ is the spin factor, $B$ is the vector of $b$ parameters,
@ -2183,7 +2183,7 @@ qmckl_exit_code qmckl_provide_jastrow_factor_ee(qmckl_context context)
ctx->jastrow.factor_ee = NULL;
}
}
/* Allocate array */
if (ctx->jastrow.factor_ee == NULL) {
@ -2264,24 +2264,24 @@ integer function qmckl_compute_factor_ee_doc_f(context, walk_num, elec_num, up_n
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 (bord_num < 0) then
info = QMCKL_INVALID_ARG_4
return
endif
factor_ee = 0.0d0
do nw =1, walk_num
do j = 1, elec_num
do i = 1, j - 1
@ -2289,23 +2289,23 @@ integer function qmckl_compute_factor_ee_doc_f(context, walk_num, elec_num, up_n
power_ser = 0.0d0
spin_fact = 1.0d0
ipar = 1
do p = 2, bord_num
x = x * ee_distance_rescaled(i,j,nw)
power_ser = power_ser + b_vector(p + 1) * x
end do
if(j <= up_num .or. i > up_num) then
spin_fact = 0.5d0
ipar = 2
endif
factor_ee(nw) = factor_ee(nw) + spin_fact * b_vector(1) * &
ee_distance_rescaled(i,j,nw) / &
(1.0d0 + b_vector(2) * &
ee_distance_rescaled(i,j,nw)) &
+ power_ser - asymp_jasb(ipar)
+ power_ser - asymp_jasb(ipar)
end do
end do
end do
@ -2417,8 +2417,8 @@ qmckl_exit_code qmckl_compute_factor_ee_hpc (
ipar = 1;
}
factor_ee[nw] += spin_fact * b_vector[0] *
x1 / (1.0 + b_vector[1] * x1)
factor_ee[nw] += spin_fact * b_vector[0] *
x1 / (1.0 + b_vector[1] * x1)
- asymp_jasb[ipar] + power_ser;
}
@ -2624,7 +2624,7 @@ qmckl_exit_code qmckl_provide_jastrow_factor_ee_deriv_e(qmckl_context context)
ctx->jastrow.factor_ee_deriv_e = NULL;
}
}
/* Allocate array */
if (ctx->jastrow.factor_ee_deriv_e == NULL) {
@ -2804,11 +2804,11 @@ qmckl_exit_code qmckl_compute_factor_ee_deriv_e_hpc(
if (context == QMCKL_NULL_CONTEXT) {
return QMCKL_INVALID_CONTEXT;
}
}
if (walk_num <= 0) {
return QMCKL_INVALID_ARG_2;
}
}
if (elec_num <= 0) {
return QMCKL_INVALID_ARG_3;
@ -2818,7 +2818,7 @@ qmckl_exit_code qmckl_compute_factor_ee_deriv_e_hpc(
return QMCKL_INVALID_ARG_4;
}
for (int nw = 0; nw < walk_num; ++nw) {
for (int ii = 0; ii < 4; ++ii) {
for (int j = 0; j < elec_num; ++j) {
@ -2826,7 +2826,7 @@ qmckl_exit_code qmckl_compute_factor_ee_deriv_e_hpc(
}
}
}
third = 1.0 / 3.0;
for (int nw = 0; nw < walk_num; ++nw) {
@ -2836,14 +2836,14 @@ qmckl_exit_code qmckl_compute_factor_ee_deriv_e_hpc(
if (fabs(x) < 1.0e-18) continue;
for (int ii = 0; ii < 3; ++ii){
pow_ser_g[ii] = 0.0;
}
}
spin_fact = 1.0;
den = 1.0 + b_vector[1] * x;
invden = 1.0 / den;
invden2 = invden * invden;
invden3 = invden2 * invden;
xinv = 1.0 / (x + 1.0e-18);
dx[0] = ee_distance_rescaled_deriv_e[0 \
+ j * 4 + i * 4 * elec_num \
+ nw * 4 * elec_num * elec_num];
@ -2891,7 +2891,7 @@ qmckl_exit_code qmckl_compute_factor_ee_deriv_e_hpc(
}
}
}
return QMCKL_SUCCESS;
}
#+end_src
@ -2956,7 +2956,7 @@ integer(c_int32_t) function qmckl_compute_factor_ee_deriv_e_doc &
end function qmckl_compute_factor_ee_deriv_e_doc
#+end_src
#+begin_src c :tangle (eval h_private_func) :comments org
qmckl_exit_code qmckl_compute_factor_ee_deriv_e_hpc (
const qmckl_context context,
@ -2998,9 +2998,9 @@ integer(c_int32_t) function qmckl_compute_factor_ee_deriv_e_doc &
double* const factor_ee_deriv_e ) {
#ifdef HAVE_HPC
return qmckl_compute_factor_ee_deriv_e_hpc(context, walk_num, elec_num, up_num, bord_num, b_vector, ee_distance_rescaled, ee_distance_rescaled_deriv_e, factor_ee_deriv_e );
return qmckl_compute_factor_ee_deriv_e_hpc(context, walk_num, elec_num, up_num, bord_num, b_vector, ee_distance_rescaled, ee_distance_rescaled_deriv_e, factor_ee_deriv_e );
#else
return qmckl_compute_factor_ee_deriv_e_doc(context, walk_num, elec_num, up_num, bord_num, b_vector, ee_distance_rescaled, ee_distance_rescaled_deriv_e, factor_ee_deriv_e );
return qmckl_compute_factor_ee_deriv_e_doc(context, walk_num, elec_num, up_num, bord_num, b_vector, ee_distance_rescaled, ee_distance_rescaled_deriv_e, factor_ee_deriv_e );
#endif
}
#+end_src
@ -3295,7 +3295,7 @@ integer function qmckl_compute_jastrow_asymp_jasa_f(context, aord_num, type_nucl
integer*8 :: i, j, p
double precision :: kappa_inv, x, asym_one
info = QMCKL_SUCCESS
if (context == QMCKL_NULL_CONTEXT) then
@ -3311,15 +3311,15 @@ integer function qmckl_compute_jastrow_asymp_jasa_f(context, aord_num, type_nucl
do i=1,type_nucl_num
kappa_inv = 1.0d0 / rescale_factor_en(i)
asymp_jasa(i) = a_vector(1,i) * kappa_inv / (1.0d0 + a_vector(2,i) * kappa_inv)
x = kappa_inv
do p = 2, aord_num
x = x * kappa_inv
asymp_jasa(i) = asymp_jasa(i) + a_vector(p+1, i) * x
end do
end do
end function qmckl_compute_jastrow_asymp_jasa_f
@ -3371,7 +3371,7 @@ qmckl_exit_code qmckl_compute_jastrow_asymp_jasa (
for (int i = 0 ; i <= type_nucl_num; ++i) {
const double kappa_inv = 1.0 / rescale_factor_en[i];
asymp_jasa[i] = a_vector[aord_num*i] * kappa_inv / (1.0 + a_vector[1 + aord_num*i] * kappa_inv);
double x = kappa_inv;
for (int p = 1; p < aord_num; ++p){
x *= kappa_inv;
@ -3394,7 +3394,7 @@ qmckl_exit_code qmckl_compute_jastrow_asymp_jasa (
const int64_t type_nucl_num,
const double* a_vector,
const double* rescale_factor_en,
double* const asymp_jasa );
double* const asymp_jasa );
#+end_src
*** Test
@ -3577,7 +3577,7 @@ qmckl_exit_code qmckl_provide_jastrow_factor_en(qmckl_context context)
if (rc != QMCKL_SUCCESS) {
return rc;
}
ctx->jastrow.factor_en_date = ctx->date;
}
@ -3829,7 +3829,7 @@ qmckl_exit_code qmckl_compute_factor_en (
const double* a_vector,
const double* en_distance_rescaled,
const double* asymp_jasa,
double* const factor_en );
double* const factor_en );
#+end_src
*** Test
@ -6501,7 +6501,7 @@ integer function qmckl_compute_een_rescaled_n_f( &
integer*8 , intent(in) :: type_nucl_vector(nucl_num)
integer*8 , intent(in) :: cord_num
double precision , intent(in) :: rescale_factor_en(type_nucl_num)
double precision , intent(in) :: en_distance(elec_num,nucl_num,walk_num)
double precision , intent(in) :: en_distance(nucl_num,elec_num,walk_num)
double precision , intent(out) :: een_rescaled_n(elec_num,nucl_num,0:cord_num,walk_num)
double precision :: x
integer*8 :: i, a, k, l, nw
@ -6536,16 +6536,16 @@ integer function qmckl_compute_een_rescaled_n_f( &
! Prepare table of exponentiated distances raised to appropriate power
een_rescaled_n = 0.0d0
do nw = 1, walk_num
! prepare the actual een table
een_rescaled_n(:, :, 0, nw) = 1.0d0
do a = 1, nucl_num
do i = 1, elec_num
een_rescaled_n(i, a, 1, nw) = dexp(-rescale_factor_en(type_nucl_vector(a)) * en_distance(i, a, nw))
een_rescaled_n(i, a, 1, nw) = dexp(-rescale_factor_en(type_nucl_vector(a)) * en_distance(a, i, nw))
end do
end do
do l = 2, cord_num
do a = 1, nucl_num
do i = 1, elec_num
@ -6604,7 +6604,7 @@ qmckl_exit_code qmckl_compute_een_rescaled_n (
for (int i = 0; i < elec_num; ++i) {
een_rescaled_n[i + a*elec_num + nw * elec_num*nucl_num*(cord_num+1)] = 1.0;
een_rescaled_n[i + a*elec_num + elec_num*nucl_num + nw*elec_num*nucl_num*(cord_num+1)] =
exp(-rescale_factor_en[type_nucl_vector[a]] * en_distance[i + a*elec_num + nw*elec_num*nucl_num]);
exp(-rescale_factor_en[type_nucl_vector[a]] * en_distance[a + i*nucl_num + nw*elec_num*nucl_num]);
}
}
@ -6612,7 +6612,7 @@ qmckl_exit_code qmckl_compute_een_rescaled_n (
for (int a = 0; a < nucl_num; ++a) {
for (int i = 0; i < elec_num; ++i) {
een_rescaled_n[i + a*elec_num + l*elec_num*nucl_num + nw*elec_num*nucl_num*(cord_num+1)] =
een_rescaled_n[i + a*elec_num + (l-1)*elec_num*nucl_num + nw*elec_num*nucl_num*(cord_num+1)] *
een_rescaled_n[i + a*elec_num + (l-1)*elec_num*nucl_num + nw*elec_num*nucl_num*(cord_num+1)] *
een_rescaled_n[i + a*elec_num + elec_num*nucl_num + nw*elec_num*nucl_num*(cord_num+1)];
}
}
@ -6671,7 +6671,7 @@ qmckl_exit_code qmckl_compute_een_rescaled_n (
end function qmckl_compute_een_rescaled_n
#+end_src
# #+CALL: generate_c_header(table=qmckl_factor_een_rescaled_n_args,rettyp=get_value("CRetType"),fname=get_value("Name"))
#+begin_src c :comments org :tangle (eval h_private_func) :noweb yes :exports none
@ -6917,7 +6917,7 @@ integer function qmckl_compute_factor_een_rescaled_n_deriv_e_f( &
double precision , intent(in) :: rescale_factor_en(type_nucl_num)
double precision , intent(in) :: coord_ee(elec_num,3,walk_num)
double precision , intent(in) :: coord_en(nucl_num,3)
double precision , intent(in) :: en_distance(elec_num,nucl_num,walk_num)
double precision , intent(in) :: en_distance(nucl_num,elec_num,walk_num)
double precision , intent(in) :: een_rescaled_n(elec_num,nucl_num,0:cord_num,walk_num)
double precision , intent(out) :: een_rescaled_n_deriv_e(elec_num,4,nucl_num,0:cord_num,walk_num)
double precision,dimension(:,:,:),allocatable :: elnuc_dist_deriv_e
@ -6960,7 +6960,7 @@ integer function qmckl_compute_factor_een_rescaled_n_deriv_e_f( &
! prepare the actual een table
do a = 1, nucl_num
do i = 1, elec_num
ria_inv = 1.0d0 / en_distance(i, a, nw)
ria_inv = 1.0d0 / en_distance(a, i, nw)
do ii = 1, 3
elnuc_dist_deriv_e(ii, i, a) = (coord_ee(i, ii, nw) - coord_en(a, ii)) * ria_inv
end do
@ -7441,7 +7441,7 @@ qmckl_exit_code qmckl_provide_dtmp_c(qmckl_context context)
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT;
}
qmckl_exit_code rc;
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
assert (ctx != NULL);
@ -7487,7 +7487,7 @@ qmckl_exit_code qmckl_provide_dtmp_c(qmckl_context context)
ctx->jastrow.een_rescaled_e_deriv_e,
ctx->jastrow.een_rescaled_n,
ctx->jastrow.dtmp_c);
if (rc != QMCKL_SUCCESS) {
return rc;
}

6
org/qmckl_point.org
View File

@ -19,7 +19,7 @@ coordinates of all the walkers.
#include <stdbool.h>
#include "qmckl_blas_private_type.h"
#+end_src
#+begin_src c :tangle (eval h_private_func)
#ifndef QMCKL_POINT_HPF
#define QMCKL_POINT_HPF
@ -318,8 +318,8 @@ qmckl_set_point (qmckl_context context,
if (num != ctx->point.num) {
if (ctx->point.coord.data != NULL) {
rc = qmckl_matrix_free(context, &(ctx->point.coord));
assert (rc == QMCKL_SUCCESS);
rc = qmckl_matrix_free(context, &(ctx->point.coord));
assert (rc == QMCKL_SUCCESS);
}
ctx->point.coord = qmckl_matrix_alloc(context, num, 3);