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

Fixed bug and working on tests for MOs. #32

This commit is contained in:
v1j4y 2021-09-27 18:17:21 +02:00
parent ff5e7882d0
commit a055732227

View File

@ -306,6 +306,9 @@ qmckl_exit_code qmckl_get_mo_basis_vgl(qmckl_context context, double* const mo_v
rc = qmckl_provide_ao_vgl(context);
if (rc != QMCKL_SUCCESS) return rc;
rc = qmckl_provide_mo_vgl(context);
if (rc != QMCKL_SUCCESS) return rc;
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
@ -333,11 +336,11 @@ qmckl_exit_code qmckl_get_mo_basis_vgl(qmckl_context context, double* const mo_v
*** Provide
#+begin_src c :comments org :tangle (eval h_private_func) :noweb yes :exports none
qmckl_exit_code qmckl_provide_mo_basis_vgl(qmckl_context context);
qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context);
#+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_provide_mo_basis_vgl(qmckl_context context)
qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context)
{
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
@ -439,7 +442,7 @@ integer function qmckl_compute_mo_basis_gaussian_vgl_f(context, &
double precision , intent(in) :: ao_vgl(ao_num,elec_num,walk_num,5)
double precision , intent(in) :: coef_normalized(mo_num,ao_num)
double precision , intent(out) :: mo_vgl(mo_num,elec_num,walk_num,5)
logical :: TransA, TransB
logical*8 :: TransA, TransB
double precision :: alpha, beta
integer :: info_qmckl_dgemm_value
integer :: info_qmckl_dgemm_Gx
@ -456,46 +459,61 @@ integer function qmckl_compute_mo_basis_gaussian_vgl_f(context, &
beta = 0.0d0
info = QMCKL_SUCCESS
info_qmckl_dgemm_value = QMCKL_SUCCESS
info_qmckl_dgemm_Gx = QMCKL_SUCCESS
info_qmckl_dgemm_Gy = QMCKL_SUCCESS
info_qmckl_dgemm_Gz = QMCKL_SUCCESS
info_qmckl_dgemm_lap = QMCKL_SUCCESS
! Don't compute exponentials when the result will be almost zero.
! TODO : Use numerical precision here
cutoff = -dlog(1.d-15)
M = mo_num * 1_8
N = ao_num * 1_8
K = 1_8
M = 1_8
N = 1_8
K = mo_num * 1_8
do iwalk = 1, walk_num
do ielec = 1, elec_num
! Value
TransA = .True.
TransB = .True.
info_qmckl_dgemm_value = qmckl_dgemm(context,TransA, TransB, M, N, K, alpha, &
ao_vgl(:, :, iwalk, 1), size(ao_vgl,1) * 1_8, &
ao_vgl(:, ielec, iwalk, 1), size(ao_vgl,1) * 1_8, &
coef_normalized,size(coef_normalized,1) * 1_8, &
beta, &
mo_vgl(:,:,iwalk,1),size(mo_vgl,1) * 1_8)
mo_vgl(:,ielec,iwalk,1),1_8)
! Grad_x
TransA = .True.
TransB = .True.
info_qmckl_dgemm_Gx = qmckl_dgemm(context,TransA, TransB, M, N, K, alpha, &
ao_vgl(:, :, iwalk, 2), size(ao_vgl,1) * 1_8, &
ao_vgl(:, ielec, iwalk, 2), size(ao_vgl,1) * 1_8, &
coef_normalized,size(coef_normalized,1) * 1_8, &
beta, &
mo_vgl(:,:,iwalk,2),size(mo_vgl,1) * 1_8)
mo_vgl(:,ielec,iwalk,2),1_8)
! Grad_y
TransA = .True.
TransB = .True.
info_qmckl_dgemm_Gy = qmckl_dgemm(context,TransA, TransB, M, N, K, alpha, &
ao_vgl(:, :, iwalk, 3), size(ao_vgl,1) * 1_8, &
ao_vgl(:, ielec, iwalk, 3), size(ao_vgl,1) * 1_8, &
coef_normalized,size(coef_normalized,1) * 1_8, &
beta, &
mo_vgl(:,:,iwalk,3),size(mo_vgl,1) * 1_8)
mo_vgl(:,ielec,iwalk,3),1_8)
! Grad_z
TransA = .True.
TransB = .True.
info_qmckl_dgemm_Gz = qmckl_dgemm(context,TransA, TransB, M, N, K, alpha, &
ao_vgl(:, :, iwalk, 4), size(ao_vgl,1) * 1_8, &
ao_vgl(:, ielec, iwalk, 4), size(ao_vgl,1) * 1_8, &
coef_normalized,size(coef_normalized,1) * 1_8, &
beta, &
mo_vgl(:,:,iwalk,4),size(mo_vgl,1) * 1_8)
mo_vgl(:,ielec,iwalk,4),1_8)
! Lapl_z
info_qmckl_dgemm_lap = qmckl_dgemm(context,TransA, TransB, M, N, K, alpha, &
ao_vgl(:, :, iwalk, 5), size(ao_vgl,1) * 1_8, &
TransA = .True.
TransB = .True.
info_qmckl_dgemm_lap = qmckl_dgemm(context, TransA, TransB, M, N, K, alpha, &
ao_vgl(:, ielec, iwalk, 5), size(ao_vgl,1) * 1_8, &
coef_normalized,size(coef_normalized,1) * 1_8, &
beta, &
mo_vgl(:,:,iwalk,5),size(mo_vgl,1) * 1_8)
mo_vgl(:,ielec,iwalk,5),1_8)
end do
end do
@ -657,9 +675,6 @@ assert(rc == QMCKL_SUCCESS);
assert(qmckl_nucleus_provided(context));
//const int64_t shell_num = chbrclf_shell_num;
//const int64_t prim_num = chbrclf_prim_num;
//const int64_t ao_num = chbrclf_ao_num;
const int64_t * nucleus_index = &(chbrclf_basis_nucleus_index[0]);
const int64_t * nucleus_shell_num = &(chbrclf_basis_nucleus_shell_num[0]);
const int32_t * shell_ang_mom = &(chbrclf_basis_shell_ang_mom[0]);
@ -736,9 +751,30 @@ double ao_vgl[5][walk_num][elec_num][chbrclf_ao_num];
rc = qmckl_get_ao_vgl(context, &(ao_vgl[0][0][0][0]));
assert (rc == QMCKL_SUCCESS);
rc = qmckl_get_mo_basis_vgl(context, &(ao_vgl[0][0][0][0]));
/* Set up MO data */
rc = qmckl_set_mo_basis_type(context, typ);
assert (rc == QMCKL_SUCCESS);
const int64_t mo_num = chbrclf_ao_num;
rc = qmckl_set_mo_basis_mo_num(context, mo_num);
assert (rc == QMCKL_SUCCESS);
double mo_coefficient[chbrclf_ao_num][mo_num];
const int64_t min_num = (((chbrclf_ao_num) < (mo_num)) ? (chbrclf_ao_num) : (mo_num));
for (int i=0; i<min_num; ++i) {
mo_coefficient[i][i] = 1.0;
}
rc = qmckl_set_mo_basis_coefficient(context, &(mo_coefficient[0][0]));
assert (rc == QMCKL_SUCCESS);
assert(qmckl_ao_basis_provided(context));
double mo_vgl[5][walk_num][elec_num][chbrclf_ao_num];
rc = qmckl_get_mo_basis_vgl(context, &(mo_vgl[0][0][0][0]));
assert (rc == QMCKL_SUCCESS);
//printf("\n");
//printf(" ao_vgl ao_vgl[0][0][26][219] %25.15e\n", ao_vgl[0][0][26][219]);
//printf(" ao_vgl ao_vgl[1][0][26][219] %25.15e\n", ao_vgl[1][0][26][219]);