TREX/qmckl
1
0
Fork 0
mirror of https://github.com/TREX-CoE/qmckl.git synced 2025-04-30 04:15:00 +02:00
Code Issues Projects Releases Wiki Activity

added more tests

Browse Source
This commit is contained in:
Emiel Slootman 2024-12-10 09:46:21 +01:00
parent c0f4d095d3
commit 0ee92ca5aa
1 changed files with 275 additions and 40 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

315
org/qmckl_jastrow_champ_single.org
View File

@ -4597,7 +4597,6 @@ integer function qmckl_compute_jastrow_champ_delta_p_gl_doc_f( &
double precision , intent(out) :: delta_p_gl(4,elec_num, nucl_num,0:cord_num, 0:cord_num-1, walk_num)
double precision :: delta_e_gl(4, elec_num, 0:cord_num, walk_num)
double precision :: delta_n_gl(4, nucl_num, 0:cord_num, walk_num)
double precision :: delta_e_gl_2(elec_num, 0:cord_num, walk_num)
double precision :: een_rescaled_e_gl_2(elec_num, elec_num, 0:cord_num, walk_num)
@ -4628,7 +4627,6 @@ integer function qmckl_compute_jastrow_champ_delta_p_gl_doc_f( &
if (cord_num == 0) return
delta_e_gl(:,:,:,:) = een_rescaled_single_e_gl(:,:,:,:) - een_rescaled_e_gl(num, :, :, :, :)
delta_n_gl(:,:,:,:) = een_rescaled_single_n_gl(:,:,:,:) - een_rescaled_n_gl(num, :, :, :, :)
delta_e_gl(:, num, :, :) = 0.0d0
een_rescaled_delta_n(:,:,:) = een_rescaled_single_n(:,:,:) - een_rescaled_n(num, :, :, :)
@ -4640,60 +4638,60 @@ integer function qmckl_compute_jastrow_champ_delta_p_gl_doc_f( &
! First calculate delta p
do nw=1, walk_num
do i=0, cord_num-1
do m=0, cord_num-1
do k = 1, 4
!do a = 1, nucl_num
!do l=0, cord_num
!do j = 1, elec_num
!delta_p_gl(k,num,a,l,i,nw)= delta_p_gl(k,num,a,l,i,nw) + delta_e_gl(k,j,i,nw) * een_rescaled_n(j,a,l,nw)
!delta_p_gl(k,j,a,l,i,nw) = delta_p_gl(k,j,a,l,i,nw) + delta_e_gl(k,j,i,nw) * een_rescaled_n(num,a,l,nw)
do a = 1, nucl_num
do l=0, cord_num
do j = 1, elec_num
delta_p_gl(k,num,a,l,m,nw) = delta_p_gl(k,num,a,l,m,nw) + delta_e_gl(k,j,m,nw) * een_rescaled_n(j,a,l,nw)
delta_p_gl(k,j,a,l,m,nw) = delta_p_gl(k,j,a,l,m,nw) + delta_e_gl(k,j,m,nw) * een_rescaled_n(num,a,l,nw)
!delta_p_gl(k,j,a,l,i,nw) = delta_p_gl(k,j,a,l,i,nw) + een_rescaled_e_gl(j,k, num,i,nw) * een_rescaled_delta_n(a,l,nw)
delta_p_gl(k,j,a,l,m,nw) = delta_p_gl(k,j,a,l,m,nw) + een_rescaled_e_gl(j,k, num,m,nw) * een_rescaled_delta_n(a,l,nw)
!delta_p_gl(k,j,a,l,i,nw) = delta_p_gl(k,j,a,l,i,nw) + delta_e_gl(k, j,i,nw) * een_rescaled_delta_n(a,l,nw)
delta_p_gl(k,j,a,l,m,nw) = delta_p_gl(k,j,a,l,m,nw) + delta_e_gl(k, j,m,nw) * een_rescaled_delta_n(a,l,nw)
!end do
!end do
!end do
end do
end do
end do
delta_e_gl_2(:,:,:) = delta_e_gl(k, :,:,:)
een_rescaled_e_gl_2(:,:,:,:) = een_rescaled_e_gl(:,k, :,:,:)
!delta_e_gl_2(:,:,:) = delta_e_gl(k, :,:,:)
!een_rescaled_e_gl_2(:,:,:,:) = een_rescaled_e_gl(:,k, :,:,:)
info = qmckl_dgemm(context, 'T', 'N', 1, nucl_num * (cord_num+1), elec_num, 1.0d0, &
delta_e_gl_2(1,i,nw),elec_num, &
een_rescaled_n(1,1,0,nw),elec_num, &
0.0d0, &
delta_c(1,0,i,nw),1)
!info = qmckl_dgemm(context, 'T', 'N', 1, nucl_num * (cord_num+1), elec_num, 1.0d0, &
! delta_e_gl_2(1,m,nw),elec_num, &
! een_rescaled_n(1,1,0,nw),elec_num, &
! 0.0d0, &
! delta_c(1,0,m,nw),1)
info = qmckl_dgemm(context, 'N', 'N', elec_num, nucl_num * (cord_num+1), 1, 1.0d0, &
delta_e_gl_2(1,i,nw),elec_num, &
een_rescaled_n(num,1,0,nw),elec_num, &
0.0d0, &
delta_c2(1,1,0,i,nw),elec_num)
!info = qmckl_dgemm(context, 'N', 'N', elec_num, nucl_num * (cord_num+1), 1, 1.0d0, &
! delta_e_gl_2(1,m,nw),elec_num, &
! een_rescaled_n(num,1,0,nw),elec_num, &
! 0.0d0, &
! delta_c2(1,1,0,m,nw),elec_num)
delta_c2(num,:,:,i,nw) = delta_c2(num,:,:,i,nw) + delta_c(:,:,i,nw)
delta_p_gl(k,:,:,:,i,nw) = delta_c2(:,:,:,i,nw)
!delta_c2(num,:,:,m,nw) = delta_c2(num,:,:,m,nw) + delta_c(:,:,m,nw)
!delta_p_gl(k,:,:,:,m,nw) = delta_c2(:,:,:,m,nw)
info = qmckl_dgemm(context, 'N', 'T', elec_num, nucl_num * (cord_num+1), 1, 1.0d0, &
een_rescaled_e_gl_2(:,num,i,nw),elec_num, &
een_rescaled_delta_n(:,:,nw),nucl_num* (cord_num+1), &
0.0d0, &
delta_c2(:,:,:,i,nw),elec_num)
!info = qmckl_dgemm(context, 'N', 'T', elec_num, nucl_num * (cord_num+1), 1, 1.0d0, &
! een_rescaled_e_gl_2(:,num,m,nw),elec_num, &
! een_rescaled_delta_n(:,:,nw),nucl_num* (cord_num+1), &
! 0.0d0, &
! delta_c2(:,:,:,m,nw),elec_num)
delta_p_gl(k,:,:,:,i,nw) = delta_p_gl(k,:,:,:,i,nw) + delta_c2(:,:,:,i,nw)
!delta_p_gl(k,:,:,:,m,nw) = delta_p_gl(k,:,:,:,m,nw) + delta_c2(:,:,:,m,nw)
info = qmckl_dgemm(context, 'N', 'T', elec_num, nucl_num * (cord_num+1), 1, 1.0d0, &
delta_e_gl_2(:,i,nw),elec_num, &
een_rescaled_delta_n(:,:,nw),nucl_num* (cord_num+1), &
0.0d0, &
delta_c2(:,:,:,i,nw),elec_num)
!info = qmckl_dgemm(context, 'N', 'T', elec_num, nucl_num * (cord_num+1), 1, 1.0d0, &
! delta_e_gl_2(:,m,nw),elec_num, &
! een_rescaled_delta_n(:,:,nw),nucl_num* (cord_num+1), &
! 0.0d0, &
! delta_c2(:,:,:,m,nw),elec_num)
delta_p_gl(k,:,:,:,i,nw) = delta_p_gl(k,:,:,:,i,nw) + delta_c2(:,:,:,i,nw)
!delta_p_gl(k,:,:,:,m,nw) = delta_p_gl(k,:,:,:,m,nw) + delta_c2(:,:,:,m,nw)
end do
@ -4877,7 +4875,7 @@ 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++){
printf("p_gl_new[%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("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, i, k, delta_p_gl[nw][l][m][a][i][k]);
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][i][k])) < 1.e-12);
@ -5307,6 +5305,54 @@ qmckl_compute_jastrow_champ_factor_single_een_gl (const qmckl_context context,
end function qmckl_compute_jastrow_champ_factor_single_een_gl_doc
#+end_src
** test
#+begin_src c :tangle (eval c_test)
/* Check if Jastrow is properly initialized */
assert(qmckl_jastrow_champ_provided(context));
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][elec_num][4];
rc = qmckl_get_jastrow_champ_factor_een_gl(context, &een_gl_old[0][0][0], walk_num*elec_num*4);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_get_electron_coord(context, 'N', &coords[0][0][0], walk_num*elec_num*3);
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);
assert (rc == QMCKL_SUCCESS);
coords[0][2][0] = new_coords[0];
coords[0][2][1] = new_coords[1];
coords[0][2][2] = new_coords[2];
rc = qmckl_set_point(context, 'N', elec_num, &coords[0][0][0], walk_num*elec_num*3);
assert (rc == QMCKL_SUCCESS);
double een_gl_new[walk_num][elec_num][4];
rc = qmckl_get_jastrow_champ_factor_een_gl(context, &een_gl_new[0][0][0], walk_num*elec_num*4);
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++) {
assert(fabs((een_gl_new[nw][i][m] - een_gl_new[nw][i][m]) - delta_een_gl[nw][i][m]) < 1.e-12);
}
}
}
#+end_src
* ee distance rescaled single point gl
** Get
@ -5559,6 +5605,54 @@ qmckl_exit_code qmckl_compute_ee_rescaled_single_gl (
}
#+end_src
** test
#+begin_src c :tangle (eval c_test)
/* Check if Jastrow is properly initialized */
assert(qmckl_jastrow_champ_provided(context));
rc = qmckl_set_point(context, 'N', elec_num, elec_coord, walk_num*elec_num*3);
assert(rc == QMCKL_SUCCESS);
double ee_rescaled_gl[walk_num][elec_num][elec_num][4];
rc = qmckl_get_jastrow_champ_ee_distance_rescaled_gl(context, &ee_rescaled_gl[0][0][0][0]);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_get_electron_coord(context, 'N', &coords[0][0][0], walk_num*elec_num*3);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_set_single_point(context, 'N', 2, new_coords, 3);
assert (rc == QMCKL_SUCCESS);
double single_ee_rescaled_gl[walk_num][elec_num][4];
rc = qmckl_get_ee_rescaled_single_gl(context, &single_ee_rescaled_gl[0][0][0]);
assert (rc == QMCKL_SUCCESS);
coords[0][2][0] = new_coords[0];
coords[0][2][1] = new_coords[1];
coords[0][2][2] = new_coords[2];
rc = qmckl_set_point(context, 'N', elec_num, &coords[0][0][0], walk_num*elec_num*3);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_get_jastrow_champ_ee_distance_rescaled_gl(context, &ee_rescaled_gl[0][0][0][0]);
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++) {
assert(fabs(ee_rescaled_gl[nw][2][i][m] - single_ee_rescaled_gl[nw][i][m]) < 1.e-12);
}
}
}
#+end_src
* en distance rescaled single point gl
** Get
@ -5859,6 +5953,53 @@ qmckl_exit_code qmckl_compute_en_rescaled_single_gl (
end function qmckl_compute_en_rescaled_single_gl_doc
#+end_src
** test
#+begin_src c :tangle (eval c_test)
/* Check if Jastrow is properly initialized */
assert(qmckl_jastrow_champ_provided(context));
rc = qmckl_set_point(context, 'N', elec_num, elec_coord, walk_num*elec_num*3);
assert(rc == QMCKL_SUCCESS);
double en_rescaled_gl[walk_num][nucl_num][elec_num][4];
rc = qmckl_get_electron_en_distance_rescaled_gl(context, &en_rescaled_gl[0][0][0][0]);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_get_electron_coord(context, 'N', &coords[0][0][0], walk_num*elec_num*3);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_set_single_point(context, 'N', 2, new_coords, 3);
assert (rc == QMCKL_SUCCESS);
double single_en_rescaled_gl[walk_num][nucl_num][4];
rc = qmckl_get_en_rescaled_single_gl(context, &single_en_rescaled_gl[0][0][0]);
assert (rc == QMCKL_SUCCESS);
coords[0][2][0] = new_coords[0];
coords[0][2][1] = new_coords[1];
coords[0][2][2] = new_coords[2];
rc = qmckl_set_point(context, 'N', elec_num, &coords[0][0][0], walk_num*elec_num*3);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_get_electron_en_distance_rescaled_gl(context, &en_rescaled_gl[0][0][0][0]);
assert (rc == QMCKL_SUCCESS);
for (int nw = 0; nw < walk_num; nw++) {
for (int a = 0; a < nucl_num; a++) {
for (int m = 0; m < 4; m++) {
assert(fabs(en_rescaled_gl[nw][a][2][m] - single_en_rescaled_gl[nw][a][m]) < 1.e-12);
}
}
}
#+end_src
* Delta ee gl
** Get
@ -6245,6 +6386,53 @@ qmckl_compute_jastrow_champ_single_ee_gl (const qmckl_context context,
}
#+end_src
** test
#+begin_src c :tangle (eval c_test)
/* Check if Jastrow is properly initialized */
assert(qmckl_jastrow_champ_provided(context));
rc = qmckl_set_point(context, 'N', elec_num, elec_coord, walk_num*elec_num*3);
assert(rc == QMCKL_SUCCESS);
double ee_gl_old[walk_num][elec_num][4];
rc = qmckl_get_jastrow_champ_factor_ee_gl(context, &ee_gl_old[0][0][0], walk_num*elec_num*4);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_get_electron_coord(context, 'N', &coords[0][0][0], walk_num*elec_num*3);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_set_single_point(context, 'N', 2, new_coords, 3);
assert (rc == QMCKL_SUCCESS);
double delta_ee_gl[walk_num][elec_num][4];
rc = qmckl_get_jastrow_champ_single_ee_gl(context, &delta_ee_gl[0][0][0], walk_num*elec_num*4);
assert (rc == QMCKL_SUCCESS);
coords[0][2][0] = new_coords[0];
coords[0][2][1] = new_coords[1];
coords[0][2][2] = new_coords[2];
rc = qmckl_set_point(context, 'N', elec_num, &coords[0][0][0], walk_num*elec_num*3);
assert (rc == QMCKL_SUCCESS);
double ee_gl_new[walk_num][elec_num][4];
rc = qmckl_get_jastrow_champ_factor_ee_gl(context, &ee_gl_new[0][0][0], walk_num*elec_num*4);
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++) {
assert(fabs((ee_gl_new[nw][i][m] - ee_gl_new[nw][i][m]) - delta_ee_gl[nw][i][m]) < 1.e-12);
}
}
}
#+end_src
* Delta en gl
** Get
@ -6617,6 +6805,53 @@ qmckl_compute_jastrow_champ_single_en_gl (const qmckl_context context,
}
#+end_src
** test
#+begin_src c :tangle (eval c_test)
/* Check if Jastrow is properly initialized */
assert(qmckl_jastrow_champ_provided(context));
rc = qmckl_set_point(context, 'N', elec_num, elec_coord, walk_num*elec_num*3);
assert(rc == QMCKL_SUCCESS);
double en_gl_old[walk_num][elec_num][4];
rc = qmckl_get_jastrow_champ_factor_en_gl(context, &en_gl_old[0][0][0], walk_num*elec_num*4);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_get_electron_coord(context, 'N', &coords[0][0][0], walk_num*elec_num*3);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_set_single_point(context, 'N', 2, new_coords, 3);
assert (rc == QMCKL_SUCCESS);
double delta_en_gl[walk_num][elec_num][4];
rc = qmckl_get_jastrow_champ_single_en_gl(context, &delta_en_gl[0][0][0], walk_num*elec_num*4);
assert (rc == QMCKL_SUCCESS);
coords[0][2][0] = new_coords[0];
coords[0][2][1] = new_coords[1];
coords[0][2][2] = new_coords[2];
rc = qmckl_set_point(context, 'N', elec_num, &coords[0][0][0], walk_num*elec_num*3);
assert (rc == QMCKL_SUCCESS);
double en_gl_new[walk_num][elec_num][4];
rc = qmckl_get_jastrow_champ_factor_en_gl(context, &en_gl_new[0][0][0], walk_num*elec_num*4);
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++) {
assert(fabs((en_gl_new[nw][i][m] - en_gl_new[nw][i][m]) - delta_en_gl[nw][i][m]) < 1.e-12);
}
}
}
#+end_src
* Accept single electron move
#+begin_src c :comments org :tangle (eval h_func) :noweb yes