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

Debug flags OK

This commit is contained in:
Anthony Scemama 2024-12-12 10:23:58 +01:00
parent 563a9e1c12
commit 77c4eb9702
2 changed files with 61 additions and 26 deletions

View File

@ -1953,7 +1953,7 @@ printf("OK\n");
#+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code
qmckl_get_forces_tmp_dc(qmckl_context context,
qmckl_get_forces_dtmp_c(qmckl_context context,
double* const forces_dtmp_c,
const int64_t size_max);
#+end_src
@ -2145,7 +2145,8 @@ integer(qmckl_exit_code) function qmckl_compute_forces_dtmp_c( &
do i = 1, elec_num
forces_dtmp_c(i,kk,a,k,m,l,nw) = 0.0
do j = 1, elec_num
forces_dtmp_c(i,kk,a,k,m,l,nw) = forces_dtmp_c(i,kk,a,k,m,l,nw) - een_rescaled_e_gl(i,kk,j,l,nw) * een_rescaled_n_gl(j,k,a,m,nw)
forces_dtmp_c(i,kk,a,k,m,l,nw) = forces_dtmp_c(i,kk,a,k,m,l,nw) - &
een_rescaled_e_gl(i,kk,j,l,nw) * een_rescaled_n_gl(j,k,a,m,nw)
end do
end do
end do
@ -2819,20 +2820,23 @@ integer(qmckl_exit_code) function qmckl_compute_forces_een_rescaled_n_gl( &
forces_een_n(i,m,a,n,l,nw) = 0.0d0
else
if (m == n) then
forces_een_n(i,m,a,n,l,nw) = forces_een_n(i,m,a,n,l,nw) + kappa_l*een_rescaled_n(i,a,l,nw)/en_distance(a,i,nw)
forces_een_n(i,m,a,n,l,nw) = forces_een_n(i,m,a,n,l,nw) + &
kappa_l*een_rescaled_n(i,a,l,nw)/en_distance(a,i,nw)
end if
if (m < 4) then
forces_een_n(i,m,a,n,l,nw) = forces_een_n(i,m,a,n,l,nw) - &
een_rescaled_n_gl(i,m,a,l,nw) * een_rescaled_n_gl(i,n,a,l,nw) /(kappa_l*een_rescaled_n(i,a,l,nw)*en_distance(a,i,nw)) - &
een_rescaled_n_gl(i,m,a,l,nw) * een_rescaled_n_gl(i,n,a,l,nw)/een_rescaled_n(i,a,l,nw)
een_rescaled_n_gl(i,m,a,l,nw) * een_rescaled_n_gl(i,n,a,l,nw) / &
(kappa_l*een_rescaled_n(i,a,l,nw)*en_distance(a,i,nw)) - &
een_rescaled_n_gl(i,m,a,l,nw) * &
een_rescaled_n_gl(i,n,a,l,nw)/een_rescaled_n(i,a,l,nw)
else
forces_een_n(i,m,a,n,l,nw) = forces_een_n(i,m,a,n,l,nw) + &
2.0d0 * een_rescaled_n_gl(i,n,a,l,nw) /(en_distance(a,i,nw)*en_distance(a,i,nw)) - &
een_rescaled_n_gl(i,m,a,l,nw) * een_rescaled_n_gl(i,n,a,l,nw)/een_rescaled_n(i,a,l,nw)
2.0d0 * een_rescaled_n_gl(i,n,a,l,nw) / &
(en_distance(a,i,nw)*en_distance(a,i,nw)) - &
een_rescaled_n_gl(i,m,a,l,nw) * &
een_rescaled_n_gl(i,n,a,l,nw)/een_rescaled_n(i,a,l,nw)
end if
end if
end do
end do
@ -3266,14 +3270,14 @@ integer(qmckl_exit_code) function qmckl_compute_forces_jastrow_een_g( &
do i = 1, elec_num
do jj = 1, 3
forces_jastrow_een_g(i, ii, a, jj, nw) = forces_jastrow_een_g(i, ii, a, jj, nw) + (&
tmp_c(i,a, m, k,nw) * forces_een_n(i,ii,a,jj ,m+l,nw) + &
forces_tmp_c(i,jj,a, m, k,nw) * een_rescaled_n_gl(i,ii,a,m+l,nw) + &
tmp_c(i,a, m+l,k,nw) * forces_een_n(i,ii,a,jj, m, nw) + &
forces_tmp_c(i,jj,a, m+l,k,nw) * een_rescaled_n_gl(i,ii,a,m, nw) + &
dtmp_c(i,ii,a, m, k,nw) * -een_rescaled_n_gl(i,jj,a,m+l,nw) + &
forces_dtmp_c(i,ii,a,jj,m, k,nw) * een_rescaled_n(i,a, m+l,nw) + &
dtmp_c(i,ii,a, m+l,k,nw) * -een_rescaled_n_gl(i,jj,a,m, nw) + &
forces_dtmp_c(i,ii,a,jj,m+l,k,nw) * een_rescaled_n(i,a, m, nw) &
tmp_c(i,a, m, k,nw) * forces_een_n(i,ii,a,jj ,m+l,nw) &
+ forces_tmp_c(i,jj,a, m, k,nw) * een_rescaled_n_gl(i,ii,a,m+l,nw) &
+ tmp_c(i,a, m+l,k,nw) * forces_een_n(i,ii,a,jj, m, nw) &
+ forces_tmp_c(i,jj,a, m+l,k,nw) * een_rescaled_n_gl(i,ii,a,m, nw) &
- dtmp_c(i,ii,a, m, k,nw) * een_rescaled_n_gl(i,jj,a,m+l,nw) &
+ forces_dtmp_c(i,ii,a,jj,m, k,nw) * een_rescaled_n(i,a, m+l,nw) &
- dtmp_c(i,ii,a, m+l,k,nw) * een_rescaled_n_gl(i,jj,a,m, nw) &
+ forces_dtmp_c(i,ii,a,jj,m+l,k,nw) * een_rescaled_n(i,a, m, nw) &
) * cn
end do
@ -3637,15 +3641,15 @@ integer(qmckl_exit_code) function qmckl_compute_forces_jastrow_een_l( &
do ii = 1, 3
do i = 1, elec_num
forces_jastrow_een_l(a, ii, nw) = forces_jastrow_een_l(a, ii, nw) + (&
tmp_c(i,a, m, k,nw) * forces_een_n(i,4,a,ii ,m+l,nw) + &
forces_tmp_c(i,ii,a, m, k,nw) * een_rescaled_n_gl(i,4,a,m+l,nw) + &
tmp_c(i,a, m+l,k,nw) * forces_een_n(i,4,a,ii, m, nw) + &
forces_tmp_c(i,ii,a, m+l,k,nw) * een_rescaled_n_gl(i,4,a,m, nw) + &
dtmp_c(i,4,a, m, k,nw) * -een_rescaled_n_gl(i,ii,a,m+l,nw) + &
forces_dtmp_c(i,4,a,ii,m, k,nw) * een_rescaled_n(i,a, m+l,nw) + &
dtmp_c(i,4,a, m+l,k,nw) * -een_rescaled_n_gl(i,ii,a,m, nw) + &
forces_dtmp_c(i,4,a,ii,m+l,k,nw) * een_rescaled_n(i,a, m, nw) &
) * cn
tmp_c(i,a, m, k,nw) * forces_een_n(i,4,a,ii ,m+l,nw) &
+ forces_tmp_c(i,ii,a, m, k,nw) * een_rescaled_n_gl(i,4,a,m+l,nw) &
+ tmp_c(i,a, m+l,k,nw) * forces_een_n(i,4,a,ii, m, nw) &
+ forces_tmp_c(i,ii,a, m+l,k,nw) * een_rescaled_n_gl(i,4,a,m, nw) &
- dtmp_c(i,4,a, m, k,nw) * een_rescaled_n_gl(i,ii,a,m+l,nw) &
+ forces_dtmp_c(i,4,a,ii,m, k,nw) * een_rescaled_n(i,a, m+l,nw) &
- dtmp_c(i,4,a, m+l,k,nw) * een_rescaled_n_gl(i,ii,a,m, nw) &
+ forces_dtmp_c(i,4,a,ii,m+l,k,nw) * een_rescaled_n(i,a, m, nw) &
) * cn
end do
end do

View File

@ -62,6 +62,37 @@
The eN and eeN parameters are the same of all identical nuclei.
Warning: The types of nuclei use zero-based indexing.
** Reformulation of the three-body part
To accelerate the computation of the three-body part, the Jastrow
factor is re-expressed as follows:
\begin{eqnarray*}
& & \sum_{\alpha=1}^{N_{\text{nucl}}}
\sum_{i=1}^{N_{\text{elec}}}
\sum_{j=1}^{i-1}
c_{lkp\alpha} \left[ g_\text{e}({r}_{ij}) \right]^k
\left[ \left[ g_\alpha({R}_{i\alpha}) \right]^l + \left[ g_\alpha({R}_{j\alpha}) \right]^l \right]
\left[ g_\alpha({R}_{i\,\alpha}) \,
g_\alpha({R}_{j\alpha}) \right]^{m} \\
& = & \frac{1}{2}
\sum_{\alpha=1}^{N_{\text{nucl}}}
\sum_{i=1}^{N_{\text{elec}}}
\sum_{j=1}^{N_{\text{elec}}}
c_{lkp\alpha} \left[ g_\text{e}({r}_{ij}) \right]^k
\left[ \left[ g_\alpha({R}_{i\alpha}) \right]^l + \left[ g_\alpha({R}_{j\alpha}) \right]^l \right]
\left[ g_\alpha({R}_{i\,\alpha}) \, g_\alpha({R}_{j\alpha}) \right]^{m} \\
\end{eqnarray*}
& = & \frac{1}{2}
\sum_{\alpha=1}^{N_{\text{nucl}}}
\sum_{i=1}^{N_{\text{elec}}}
\sum_{j=1}^{N_{\text{elec}}}
c_{lkp\alpha} \left[ g_\text{e}({r}_{ij}) \right]^k
\left[ \left[ g_\alpha({R}_{i\alpha}) \right]^{l+m} \left[ g_\alpha({R}_{j\alpha}) \right]^{m} +
\left[ g_\alpha({R}_{i\alpha}) \right]^{l} \left[
g_\alpha({R}_{j\alpha}) \right]^{l+m} \right] \\
* Headers :noexport:
#+begin_src elisp :noexport :results none
(org-babel-lob-ingest "../tools/lib.org")