1
0
mirror of https://github.com/TREX-CoE/qmckl.git synced 2025-04-30 04:15:00 +02:00

Optimized delta_p_gl

This commit is contained in:
Anthony Scemama 2024-12-11 18:11:01 +01:00
parent 9e03e49beb
commit fde224580d

View File

@ -4227,10 +4227,10 @@ qmckl_exit_code qmckl_provide_jastrow_champ_delta_p_gl(qmckl_context context)
| ~een_rescaled_e_gl~ | ~double[walk_num][0:cord_num][elec_num][4][elec_num]~ | in | Electron-electron rescaled distances |
| ~een_rescaled_single_n_gl~ | ~double[walk_num][0:cord_num][nucl_num][4]~ | in | Electron-nucleus single rescaled distances |
| ~een_rescaled_single_e_gl~ | ~double[walk_num][0:cord_num][elec_num][4]~ | in | Electron-electron single rescaled distances |
| ~delta_p_gl~ | ~double[walk_num][0:cord_num-1][0:cord_num][nucl_num][4][elec_num]~ | out | Electron-nucleus jastrow |
| ~delta_p_gl~ | ~double[walk_num][0:cord_num-1][0:cord_num][4][nucl_num][elec_num]~ | out | Electron-nucleus jastrow |
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
integer function qmckl_compute_jastrow_champ_delta_p_gl_doc( &
integer(qmckl_exit_code) function qmckl_compute_jastrow_champ_delta_p_gl_doc( &
context, num_in, walk_num, elec_num, nucl_num, cord_num, &
een_rescaled_n, een_rescaled_e, een_rescaled_single_n, een_rescaled_single_e, &
een_rescaled_n_gl, een_rescaled_e_gl, een_rescaled_single_n_gl, een_rescaled_single_e_gl, delta_p_gl) &
@ -4248,19 +4248,17 @@ integer function qmckl_compute_jastrow_champ_delta_p_gl_doc( &
real(c_double) , intent(in) :: een_rescaled_e_gl(elec_num, 4, elec_num, 0:cord_num, walk_num)
real(c_double) , intent(in) :: een_rescaled_single_n_gl(4, nucl_num, 0:cord_num, walk_num)
real(c_double) , intent(in) :: een_rescaled_single_e_gl(4,elec_num, 0:cord_num, walk_num)
real(c_double) , intent(out) :: delta_p_gl(elec_num,4,nucl_num,0:cord_num, 0:cord_num-1, walk_num)
real(c_double) , intent(out) :: delta_p_gl(elec_num,nucl_num,4,0:cord_num, 0:cord_num-1, walk_num)
double precision :: delta_e_gl(4,elec_num)
double precision :: delta_e_gl(elec_num,4)
double precision :: delta_e_gl_2(elec_num)
double precision :: een_rescaled_e_gl_2(elec_num)
double precision :: een_rescaled_delta_n(nucl_num, 0:cord_num)
double precision :: delta_c(nucl_num,0:cord_num)
double precision :: delta_c2(elec_num, nucl_num,0:cord_num)
double precision :: een_rescaled_delta_n, fk(4)
integer*8 :: i, a, j, l, k, p, m, n, nw, num
double precision :: accu, accu2, cn
integer*8 :: LDA, LDB, LDC
double precision :: tmp
integer*8 :: LDA, LDB, LDC
num = num_in + 1
@ -4274,150 +4272,40 @@ integer function qmckl_compute_jastrow_champ_delta_p_gl_doc( &
if (info /= QMCKL_SUCCESS) return
if (cord_num == 0) return
if (cord_num == 0) then
delta_p_gl = 0.d0
return
endif
fk(1:3) = -1.d0
fk(4) = 1.d0
do nw=1, walk_num
do m=0, cord_num-1
do k = 1, 4
do j = 1, elec_num
delta_e_gl(j,k) = een_rescaled_single_e_gl(k,j,m,nw) - een_rescaled_e_gl(num, k, j, m, nw)
end do
delta_e_gl(num,k) = 0.0d0
end do
delta_c = 0.0d0
delta_c2 = 0.0d0
delta_p_gl = 0.0d0
if (.false.) then
! Compute using loops
do nw=1, walk_num
een_rescaled_delta_n(:,:) = een_rescaled_single_n(:,:,nw) - een_rescaled_n(num, :, :, nw)
do m=0, cord_num-1
delta_e_gl(:,:) = een_rescaled_single_e_gl(:,:,m,nw) - een_rescaled_e_gl(num, :, :, m, nw)
delta_e_gl(:, num) = 0.0d0
do a = 1, nucl_num
do l=0, cord_num
do l=0, cord_num
do k = 1, 4
do a = 1, nucl_num
een_rescaled_delta_n = een_rescaled_single_n(a,l,nw) - een_rescaled_n(num, a, l, nw)
tmp = fk(k) * (een_rescaled_n(num,a,l,nw) + een_rescaled_delta_n)
do j = 1, elec_num
do k = 1, 3
delta_p_gl(num,k,a,l,m,nw) = delta_p_gl(num,k,a,l,m,nw) + &
delta_e_gl(k,j) * een_rescaled_n(j,a,l,nw)
delta_p_gl(j,k,a,l,m,nw) = delta_p_gl(j,k,a,l,m,nw) - &
delta_e_gl(k,j) * een_rescaled_n(num,a,l,nw)
delta_p_gl(j,k,a,l,m,nw) = delta_p_gl(j,k,a,l,m,nw) + &
een_rescaled_e_gl(j,k,num,m,nw) * een_rescaled_delta_n(a,l)
delta_p_gl(j,k,a,l,m,nw) = delta_p_gl(j,k,a,l,m,nw) - &
delta_e_gl(k,j) * een_rescaled_delta_n(a,l)
end do
delta_p_gl(num,4,a,l,m,nw) = delta_p_gl(num,4,a,l,m,nw) + &
delta_e_gl(4,j) * een_rescaled_n(j,a,l,nw)
delta_p_gl(j,4,a,l,m,nw) = delta_p_gl(j,4,a,l,m,nw) + &
delta_e_gl(4,j) * een_rescaled_n(num,a,l,nw)
delta_p_gl(j,4,a,l,m,nw) = delta_p_gl(j,4,a,l,m,nw) + &
een_rescaled_e_gl(num,4,j,m,nw) * een_rescaled_delta_n(a,l)
delta_p_gl(j,4,a,l,m,nw) = delta_p_gl(j,4,a,l,m,nw) + &
delta_e_gl(4,j) * een_rescaled_delta_n(a,l)
delta_p_gl(j,a,k,l,m,nw) = delta_e_gl(j,k) * tmp + &
een_rescaled_e_gl(j,k,num,m,nw) * een_rescaled_delta_n
end do
do j = 1, elec_num
delta_p_gl(num,a,k,l,m,nw) = delta_p_gl(num,a,k,l,m,nw) + delta_e_gl(j,k) * een_rescaled_n(j,a,l,nw)
end do
end do
end do
end do
end do
else
! Use DGEMM
do nw=1, walk_num
een_rescaled_delta_n(:,:) = een_rescaled_single_n(:,:,nw) - een_rescaled_n(num, :, :, nw)
do m=0, cord_num-1
delta_e_gl(:,:) = een_rescaled_single_e_gl(:,:,m,nw) - een_rescaled_e_gl(num, :, :, m, nw)
delta_e_gl(:, num) = 0.0d0
do k = 1, 3
delta_e_gl_2(:) = delta_e_gl(k, :)
info = qmckl_dgemm(context, 'T', 'N', 1_8, nucl_num * (cord_num+1), elec_num, 1.0d0, &
delta_e_gl_2,elec_num, &
een_rescaled_n(1,1,0,nw),elec_num, &
0.0d0, &
delta_c,1_8)
info = qmckl_dgemm(context, 'N', 'N', elec_num, nucl_num * (cord_num+1), 1_8, -1.0d0, &
delta_e_gl_2,elec_num, &
een_rescaled_n(num,1,0,nw),elec_num, &
0.0d0, &
delta_c2(1,1,0),elec_num)
delta_c2(num,:,:) = delta_c2(num,:,:) + delta_c(:,:)
delta_p_gl(:,k,:,:,m,nw) = delta_c2(:,:,:)
info = qmckl_dgemm(context, 'N', 'T', elec_num, nucl_num * (cord_num+1), 1_8, 1.0d0, &
een_rescaled_e_gl(1,k, num,m,nw),elec_num, &
een_rescaled_delta_n,nucl_num* (cord_num+1), &
0.0d0, &
delta_c2,elec_num)
delta_p_gl(:,k,:,:,m,nw) = delta_p_gl(:,k,:,:,m,nw) + delta_c2(:,:,:)
info = qmckl_dgemm(context, 'N', 'T', elec_num, nucl_num * (cord_num+1), 1_8, -1.0d0, &
delta_e_gl_2,elec_num, &
een_rescaled_delta_n,nucl_num* (cord_num+1), &
0.0d0, &
delta_c2(:,:,:),elec_num)
delta_p_gl(:,k,:,:,m,nw) = delta_p_gl(:,k,:,:,m,nw) + delta_c2(:,:,:)
end do
k = 4
delta_e_gl_2(:) = delta_e_gl(k, :)
info = qmckl_dgemm(context, 'T', 'N', 1_8, nucl_num * (cord_num+1), elec_num, 1.0d0, &
delta_e_gl_2,elec_num, &
een_rescaled_n(1,1,0,nw),elec_num, &
0.0d0, &
delta_c,1_8)
info = qmckl_dgemm(context, 'N', 'N', elec_num, nucl_num * (cord_num+1), 1_8, 1.0d0, &
delta_e_gl_2,elec_num, &
een_rescaled_n(num,1,0,nw),elec_num, &
0.0d0, &
delta_c2,elec_num)
delta_c2(num,:,:) = delta_c2(num,:,:) + delta_c(:,:)
delta_p_gl(:,k,:,:,m,nw) = delta_c2(:,:,:)
info = qmckl_dgemm(context, 'N', 'T', elec_num, nucl_num * (cord_num+1), 1_8, 1.0d0, &
een_rescaled_e_gl(:,k, num,m,nw),elec_num, &
een_rescaled_delta_n,nucl_num* (cord_num+1), &
0.0d0, &
delta_c2,elec_num)
delta_p_gl(:,k,:,:,m,nw) = delta_p_gl(:,k,:,:,m,nw) + delta_c2(:,:,:)
info = qmckl_dgemm(context, 'N', 'T', elec_num, nucl_num * (cord_num+1), 1_8, 1.0d0, &
delta_e_gl_2,elec_num, &
een_rescaled_delta_n,nucl_num* (cord_num+1), &
0.0d0, &
delta_c2,elec_num)
delta_p_gl(:,k,:,:,m,nw) = delta_p_gl(:,k,:,:,m,nw) + delta_c2(:,:,:)
end do
end do
end if
end do
end function qmckl_compute_jastrow_champ_delta_p_gl_doc
#+end_src
@ -4512,7 +4400,7 @@ assert (rc == QMCKL_SUCCESS);
rc = qmckl_set_single_point(context, 'N', 2, new_coords, 3);
assert (rc == QMCKL_SUCCESS);
double delta_p_gl[walk_num][cord_num][cord_num+1][nucl_num][4][elec_num];
double delta_p_gl[walk_num][cord_num][cord_num+1][4][nucl_num][elec_num];
rc = qmckl_get_jastrow_champ_delta_p_gl(context, &delta_p_gl[0][0][0][0][0][0], 4*walk_num*cord_num*(cord_num+1)*nucl_num*elec_num);
assert (rc == QMCKL_SUCCESS);
@ -4534,12 +4422,12 @@ for (int nw = 0; nw < walk_num; nw++){
for (int a = 0; a < nucl_num; a++) {
for (int i = 0; i < elec_num; i++){
for (int k = 0; k < 4; k++){
if (fabs(((p_gl_new[nw][l][m][a][k][i]-p_gl_old[nw][l][m][a][k][i])-delta_p_gl[nw][l][m][a][k][i])) > 1.e-12) {
if (fabs(((p_gl_new[nw][l][m][a][k][i]-p_gl_old[nw][l][m][a][k][i])-delta_p_gl[nw][l][m][k][a][i])) > 1.e-12) {
printf("p_gl[%d][%d][%d][%d][%d][%d] = %f\n", nw, l, m, a, k, i, p_gl_new[nw][l][m][a][k][i] - p_gl_old[nw][l][m][a][k][i]);
printf("delta_p_gl[%d][%d][%d][%d][%d][%d] = %f\n", nw, l, m, a, k, i, delta_p_gl[nw][l][m][a][k][i]);
printf("delta_p_gl[%d][%d][%d][%d][%d][%d] = %f\n", nw, l, m, a, k, i, delta_p_gl[nw][l][m][k][a][i]);
}
assert(fabs(((p_gl_new[nw][l][m][a][k][i]-p_gl_old[nw][l][m][a][k][i])-delta_p_gl[nw][l][m][a][k][i])) < 1.e-12);
assert(fabs(((p_gl_new[nw][l][m][a][k][i]-p_gl_old[nw][l][m][a][k][i])-delta_p_gl[nw][l][m][k][a][i])) < 1.e-12);
}
}
}
@ -4717,7 +4605,7 @@ qmckl_exit_code qmckl_provide_jastrow_champ_single_een_gl(qmckl_context context)
| ~tmp_c~ | ~double[walk_num][0:cord_num-1][0:cord_num][nucl_num][elec_num]~ | in | vector of non-zero coefficients |
| ~dtmp_c~ | ~double[walk_num][0:cord_num-1][0:cord_num][nucl_num][4][elec_num]~ | in | vector of non-zero coefficients |
| ~delta_p~ | ~double[walk_num][0:cord_num-1][0:cord_num][nucl_num][elec_num]~ | in | vector of non-zero coefficients |
| ~delta_p_gl~ | ~double[walk_num][0:cord_num-1][0:cord_num][nucl_num][4][elec_num]~ | in | vector of non-zero coefficients |
| ~delta_p_gl~ | ~double[walk_num][0:cord_num-1][0:cord_num][4][nucl_num][elec_num]~ | in | vector of non-zero coefficients |
| ~een_rescaled_n~ | ~double[walk_num][0:cord_num][nucl_num][elec_num]~ | in | Electron-nucleus rescaled distances |
| ~een_rescaled_single_n~ | ~double[walk_num][0:cord_num][nucl_num]~ | in | Electron-nucleus single rescaled distances |
| ~een_rescaled_n_gl~ | ~double[walk_num][0:cord_num][nucl_num][4][elec_num]~ | in | Electron-nucleus rescaled distances |
@ -4725,27 +4613,28 @@ qmckl_exit_code qmckl_provide_jastrow_champ_single_een_gl(qmckl_context context)
| ~delta_een_gl~ | ~double[walk_num][elec_num][4]~ | out | Electron-nucleus jastrow |
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
integer function qmckl_compute_jastrow_champ_factor_single_een_gl_doc_f( &
integer(qmckl_exit_code) function qmckl_compute_jastrow_champ_factor_single_een_gl_doc( &
context, num_in, walk_num, elec_num, nucl_num, cord_num, &
dim_c_vector, c_vector_full, lkpm_combined_index, &
tmp_c, dtmp_c, delta_p, delta_p_gl, een_rescaled_n, een_rescaled_single_n, &
een_rescaled_n_gl, een_rescaled_single_n_gl, delta_een_gl) &
result(info)
result(info) bind(C)
use, intrinsic :: iso_c_binding
use qmckl
implicit none
integer(qmckl_context), intent(in) :: context
integer*8 , intent(in) :: num_in, walk_num, elec_num, cord_num, nucl_num, dim_c_vector
integer*8 , intent(in) :: lkpm_combined_index(dim_c_vector,4)
double precision , intent(in) :: c_vector_full(nucl_num, dim_c_vector)
double precision , intent(in) :: tmp_c(elec_num, nucl_num,0:cord_num, 0:cord_num-1, walk_num)
double precision , intent(in) :: dtmp_c(elec_num, 4, nucl_num,0:cord_num, 0:cord_num-1, walk_num)
double precision , intent(in) :: delta_p(elec_num, nucl_num,0:cord_num, 0:cord_num-1, walk_num)
double precision , intent(in) :: delta_p_gl(elec_num, 4, nucl_num,0:cord_num, 0:cord_num-1, walk_num)
double precision , intent(in) :: een_rescaled_n(elec_num, nucl_num, 0:cord_num, walk_num)
double precision , intent(in) :: een_rescaled_single_n(nucl_num, 0:cord_num, walk_num)
double precision , intent(in) :: een_rescaled_n_gl(elec_num, 4, nucl_num, 0:cord_num, walk_num)
double precision , intent(in) :: een_rescaled_single_n_gl(4, nucl_num, 0:cord_num, walk_num)
double precision , intent(out) :: delta_een_gl(4, elec_num, walk_num)
integer(c_int64_t) , intent(in), value :: num_in, walk_num, elec_num, cord_num, nucl_num, dim_c_vector
integer(c_int64_t) , intent(in) :: lkpm_combined_index(dim_c_vector,4)
real(c_double) , intent(in) :: c_vector_full(nucl_num, dim_c_vector)
real(c_double) , intent(in) :: tmp_c(elec_num, nucl_num,0:cord_num, 0:cord_num-1, walk_num)
real(c_double) , intent(in) :: dtmp_c(elec_num, 4, nucl_num,0:cord_num, 0:cord_num-1, walk_num)
real(c_double) , intent(in) :: delta_p(elec_num, nucl_num,0:cord_num, 0:cord_num-1, walk_num)
real(c_double) , intent(in) :: delta_p_gl(elec_num, nucl_num, 4, 0:cord_num, 0:cord_num-1, walk_num)
real(c_double) , intent(in) :: een_rescaled_n(elec_num, nucl_num, 0:cord_num, walk_num)
real(c_double) , intent(in) :: een_rescaled_single_n(nucl_num, 0:cord_num, walk_num)
real(c_double) , intent(in) :: een_rescaled_n_gl(elec_num, 4, nucl_num, 0:cord_num, walk_num)
real(c_double) , intent(in) :: een_rescaled_single_n_gl(4, nucl_num, 0:cord_num, walk_num)
real(c_double) , intent(out) :: delta_een_gl(elec_num, 4, walk_num)
integer*8 :: i, a, j, l, k, p, m, n, nw, kk, num
double precision :: accu, accu2, cn
@ -4780,48 +4669,49 @@ integer function qmckl_compute_jastrow_champ_factor_single_een_gl_doc_f( &
p = lkpm_combined_index(n, 3)
m = lkpm_combined_index(n, 4)
do a = 1, nucl_num
cn = c_vector_full(a, n)
if(cn == 0.d0) cycle
do i = 1, elec_num
do kk = 1, 4
delta_een_gl(kk,i,nw) = delta_een_gl(kk,i,nw) + ( &
delta_p_gl(i,kk,a,m ,k,nw) * een_rescaled_n(i,a,m+l,nw) + &
delta_p_gl(i,kk,a,m+l,k,nw) * een_rescaled_n(i,a,m ,nw) + &
do kk = 1, 4
do a = 1, nucl_num
cn = c_vector_full(a, n)
if(cn == 0.d0) cycle
do i = 1, elec_num
delta_een_gl(i,kk,nw) = delta_een_gl(i,kk,nw) + ( &
delta_p_gl(i,a,kk,m ,k,nw) * een_rescaled_n(i,a,m+l,nw) + &
delta_p_gl(i,a,kk,m+l,k,nw) * een_rescaled_n(i,a,m ,nw) + &
delta_p(i,a,m ,k,nw) * een_rescaled_n_gl(i,kk,a,m+l,nw) + &
delta_p(i,a,m+l,k,nw) * een_rescaled_n_gl(i,kk,a,m ,nw) ) * cn
end do
end do
do kk = 1, 4
delta_een_gl(kk,num,nw) = delta_een_gl(kk,num,nw) + ( &
(dtmp_c(num,kk,a,m ,k,nw) + delta_p_gl(num,kk,a,m ,k,nw)) * een_rescaled_delta_n(a,m+l,nw) + &
(dtmp_c(num,kk,a,m+l,k,nw) + delta_p_gl(num,kk,a,m+l,k,nw)) * een_rescaled_delta_n(a,m ,nw) + &
delta_een_gl(num,kk,nw) = delta_een_gl(num,kk,nw) + ( &
(dtmp_c(num,kk,a,m ,k,nw) + delta_p_gl(num,a,kk,m ,k,nw)) * een_rescaled_delta_n(a,m+l,nw) + &
(dtmp_c(num,kk,a,m+l,k,nw) + delta_p_gl(num,a,kk,m+l,k,nw)) * een_rescaled_delta_n(a,m ,nw) + &
(tmp_c(num,a,m ,k,nw) + delta_p(num,a,m ,k,nw)) * een_rescaled_delta_n_gl(kk,a,m+l,nw) + &
(tmp_c(num,a,m+l,k,nw) + delta_p(num,a,m+l,k,nw)) * een_rescaled_delta_n_gl(kk,a,m ,nw) ) &
,* cn
(tmp_c(num,a,m+l,k,nw) + delta_p(num,a,m+l,k,nw)) * een_rescaled_delta_n_gl(kk,a,m ,nw) )* cn
end do
end do
do a = 1, nucl_num
cn = c_vector_full(a, n)
if(cn == 0.d0) cycle
cn = cn + cn
do i = 1, elec_num
delta_een_gl(4,i,nw) = delta_een_gl(4,i,nw) + ( &
delta_p_gl(i,1,a,m ,k,nw) * een_rescaled_n_gl(i,1,a,m+l,nw) + &
delta_p_gl(i,1,a,m+l,k,nw) * een_rescaled_n_gl(i,1,a,m ,nw) + &
delta_p_gl(i,2,a,m ,k,nw) * een_rescaled_n_gl(i,2,a,m+l,nw) + &
delta_p_gl(i,2,a,m+l,k,nw) * een_rescaled_n_gl(i,2,a,m ,nw) + &
delta_p_gl(i,3,a,m ,k,nw) * een_rescaled_n_gl(i,3,a,m+l,nw) + &
delta_p_gl(i,3,a,m+l,k,nw) * een_rescaled_n_gl(i,3,a,m ,nw) ) * cn
delta_een_gl(i,4,nw) = delta_een_gl(i,4,nw) + ( &
delta_p_gl(i,a,1,m ,k,nw) * een_rescaled_n_gl(i,1,a,m+l,nw) + &
delta_p_gl(i,a,1,m+l,k,nw) * een_rescaled_n_gl(i,1,a,m ,nw) + &
delta_p_gl(i,a,2,m ,k,nw) * een_rescaled_n_gl(i,2,a,m+l,nw) + &
delta_p_gl(i,a,2,m+l,k,nw) * een_rescaled_n_gl(i,2,a,m ,nw) + &
delta_p_gl(i,a,3,m ,k,nw) * een_rescaled_n_gl(i,3,a,m+l,nw) + &
delta_p_gl(i,a,3,m+l,k,nw) * een_rescaled_n_gl(i,3,a,m ,nw) ) * cn
end do
delta_een_gl(4,num,nw) = delta_een_gl(4,num,nw) + ( &
(delta_p_gl(num,1,a,m ,k,nw) + dtmp_c(num,1,a,m ,k,nw)) * een_rescaled_delta_n_gl(1,a,m+l,nw) + &
(delta_p_gl(num,1,a,m+l,k,nw) + dtmp_c(num,1,a,m+l,k,nw)) * een_rescaled_delta_n_gl(1,a,m ,nw) + &
(delta_p_gl(num,2,a,m ,k,nw) + dtmp_c(num,2,a,m ,k,nw)) * een_rescaled_delta_n_gl(2,a,m+l,nw) + &
(delta_p_gl(num,2,a,m+l,k,nw) + dtmp_c(num,2,a,m+l,k,nw)) * een_rescaled_delta_n_gl(2,a,m ,nw) + &
(delta_p_gl(num,3,a,m ,k,nw) + dtmp_c(num,3,a,m ,k,nw)) * een_rescaled_delta_n_gl(3,a,m+l,nw) + &
(delta_p_gl(num,3,a,m+l,k,nw) + dtmp_c(num,3,a,m+l,k,nw)) * een_rescaled_delta_n_gl(3,a,m ,nw) ) * cn
delta_een_gl(num,4,nw) = delta_een_gl(num,4,nw) + ( &
(delta_p_gl(num,a,1,m ,k,nw) + dtmp_c(num,1,a,m ,k,nw)) * een_rescaled_delta_n_gl(1,a,m+l,nw) + &
(delta_p_gl(num,a,1,m+l,k,nw) + dtmp_c(num,1,a,m+l,k,nw)) * een_rescaled_delta_n_gl(1,a,m ,nw) + &
(delta_p_gl(num,a,2,m ,k,nw) + dtmp_c(num,2,a,m ,k,nw)) * een_rescaled_delta_n_gl(2,a,m+l,nw) + &
(delta_p_gl(num,a,2,m+l,k,nw) + dtmp_c(num,2,a,m+l,k,nw)) * een_rescaled_delta_n_gl(2,a,m ,nw) + &
(delta_p_gl(num,a,3,m ,k,nw) + dtmp_c(num,3,a,m ,k,nw)) * een_rescaled_delta_n_gl(3,a,m+l,nw) + &
(delta_p_gl(num,a,3,m+l,k,nw) + dtmp_c(num,3,a,m+l,k,nw)) * een_rescaled_delta_n_gl(3,a,m ,nw) ) * cn
end do
end do
end do
end function qmckl_compute_jastrow_champ_factor_single_een_gl_doc_f
end function qmckl_compute_jastrow_champ_factor_single_een_gl_doc
#+end_src
# #+CALL: generate_c_header(table=qmckl_factor_een_args,rettyp=qmckl_exit_code),fname=get_value("Name"))
@ -4900,76 +4790,7 @@ qmckl_compute_jastrow_champ_factor_single_een_gl (const qmckl_context context,
}
#+end_src
#+CALL: generate_c_interface(table=qmckl_factor_single_een_gl_args,rettyp=get_value("CRetType"),fname="qmckl_compute_jastrow_champ_factor_single_een_gl_doc"))
#+RESULTS:
#+begin_src f90 :tangle (eval f) :comments org :exports none
integer(c_int32_t) function qmckl_compute_jastrow_champ_factor_single_een_gl_doc &
(context, &
num, &
walk_num, &
elec_num, &
nucl_num, &
cord_num, &
dim_c_vector, &
c_vector_full, &
lkpm_combined_index, &
tmp_c, &
dtmp_c, &
delta_p, &
delta_p_gl, &
een_rescaled_n, &
een_rescaled_single_n, &
een_rescaled_n_gl, &
een_rescaled_single_n_gl, &
delta_een_gl) &
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 :: num
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 :: cord_num
integer (c_int64_t) , intent(in) , value :: dim_c_vector
real (c_double ) , intent(in) :: c_vector_full(nucl_num,dim_c_vector)
integer (c_int64_t) , intent(in) :: lkpm_combined_index(dim_c_vector,4)
real (c_double ) , intent(in) :: tmp_c(elec_num,nucl_num,0:cord_num,0:cord_num-1,walk_num)
real (c_double ) , intent(in) :: dtmp_c(elec_num,4,nucl_num,0:cord_num,0:cord_num-1,walk_num)
real (c_double ) , intent(in) :: delta_p(elec_num,nucl_num,0:cord_num,0:cord_num-1,walk_num)
real (c_double ) , intent(in) :: delta_p_gl(elec_num,4,nucl_num,0:cord_num,0:cord_num-1,walk_num)
real (c_double ) , intent(in) :: een_rescaled_n(elec_num,nucl_num,0:cord_num,walk_num)
real (c_double ) , intent(in) :: een_rescaled_single_n(nucl_num,0:cord_num,walk_num)
real (c_double ) , intent(in) :: een_rescaled_n_gl(elec_num,4,nucl_num,0:cord_num,walk_num)
real (c_double ) , intent(in) :: een_rescaled_single_n_gl(4,nucl_num,0:cord_num,walk_num)
real (c_double ) , intent(out) :: delta_een_gl(4, elec_num, walk_num)
integer(c_int32_t), external :: qmckl_compute_jastrow_champ_factor_single_een_gl_doc_f
info = qmckl_compute_jastrow_champ_factor_single_een_gl_doc_f &
(context, &
num, &
walk_num, &
elec_num, &
nucl_num, &
cord_num, &
dim_c_vector, &
c_vector_full, &
lkpm_combined_index, &
tmp_c, &
dtmp_c, &
delta_p, &
delta_p_gl, &
een_rescaled_n, &
een_rescaled_single_n, &
een_rescaled_n_gl, &
een_rescaled_single_n_gl, &
delta_een_gl)
end function qmckl_compute_jastrow_champ_factor_single_een_gl_doc
#+end_src
** test
@ -4983,7 +4804,7 @@ rc = qmckl_set_point(context, 'N', elec_num, elec_coord, walk_num*elec_num*3);
assert(rc == QMCKL_SUCCESS);
double een_gl_old[walk_num][4][elec_num];
rc = qmckl_get_jastrow_champ_factor_een_gl(context, &een_gl_old[0][0][0], walk_num*elec_num*4);
rc = qmckl_get_jastrow_champ_factor_een_gl(context, &een_gl_old[0][0][0], walk_num*4*elec_num);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_get_electron_coord(context, 'N', &coords[0][0][0], walk_num*elec_num*3);
@ -4992,8 +4813,8 @@ assert (rc == QMCKL_SUCCESS);
rc = qmckl_set_single_point(context, 'N', 2, new_coords, 3);
assert (rc == QMCKL_SUCCESS);
double delta_een_gl[walk_num][elec_num][4];
rc = qmckl_get_jastrow_champ_single_een_gl(context, &delta_een_gl[0][0][0], walk_num*elec_num*4);
double delta_een_gl[walk_num][4][elec_num];
rc = qmckl_get_jastrow_champ_single_een_gl(context, &delta_een_gl[0][0][0], walk_num*4*elec_num);
assert (rc == QMCKL_SUCCESS);
coords[0][2][0] = new_coords[0];
@ -5004,16 +4825,16 @@ rc = qmckl_set_point(context, 'N', elec_num, &coords[0][0][0], walk_num*elec_num
assert (rc == QMCKL_SUCCESS);
double een_gl_new[walk_num][4][elec_num];
rc = qmckl_get_jastrow_champ_factor_een_gl(context, &een_gl_new[0][0][0], walk_num*elec_num*4);
rc = qmckl_get_jastrow_champ_factor_een_gl(context, &een_gl_new[0][0][0], walk_num*4*elec_num);
assert (rc == QMCKL_SUCCESS);
for (int nw = 0; nw < walk_num; nw++) {
for (int i = 0; i < elec_num; i++) {
for (int m = 0; m < 4; m++) {
for (int m = 0; m < 4; m++) {
for (int i = 0; i < elec_num; i++) {
//printf("delta_een_gl[%d][%d][%d] = %f\n", nw, i, m, delta_een_gl[nw][i][m]);
//printf("een_gl_[%d][%d][%d] = %f\n", nw, m,i, een_gl_new[nw][m][i]-een_gl_old[nw][m][i]);
assert(fabs((een_gl_new[nw][m][i]- een_gl_old[nw][m][i]) - delta_een_gl[nw][i][m]) < 1.e-12);
assert(fabs((een_gl_new[nw][m][i]- een_gl_old[nw][m][i]) - delta_een_gl[nw][m][i]) < 1.e-12);
}
}
@ -6583,7 +6404,7 @@ qmckl_get_jastrow_champ_single_accept(qmckl_context context)
ctx->jastrow_champ.een_rescaled_e_gl[nw][l][ctx->single_point.num][k][i] = ctx->jastrow_champ.een_rescaled_single_e_gl[nw][l][i][k];
for (m = 0; m < ctx->jastrow_champ.cord_num; m++){
for (a = 0; a < ctx->nucleus.num; a++){
ctx->jastrow_champ.dtmp_c[nw][i][k][a][l][m] = ctx->jastrow_champ.dtmp_c[nw][i][k][a][l][m] + ctx->single_point.delta_p_gl[nw][m][l][a][i][k]
ctx->jastrow_champ.dtmp_c[nw][i][k][a][l][m] = ctx->jastrow_champ.dtmp_c[nw][i][k][a][l][m] + ctx->single_point.delta_p_gl[nw][m][l][k][a][i]
}
}
}