mirror of
https://github.com/TREX-CoE/qmckl.git
synced 2024-12-23 04:44:03 +01:00
Created a function to provide the derivative functions along with Doc. #17
This commit is contained in:
parent
55d749f5c6
commit
cb70c1f568
@ -1105,6 +1105,360 @@ end function qmckl_distance_rescaled_f
|
||||
#+end_src
|
||||
|
||||
*** Test :noexport:
|
||||
* Rescaled Distance Derivatives
|
||||
|
||||
** ~qmckl_distance_rescaled_deriv_e~
|
||||
:PROPERTIES:
|
||||
:Name: qmckl_distance_rescaled_deriv_e
|
||||
:CRetType: qmckl_exit_code
|
||||
:FRetType: qmckl_exit_code
|
||||
:END:
|
||||
|
||||
~qmckl_distance_rescaled_deriv_e~ computes the matrix of the gradient and laplacian of the
|
||||
rescaled distance with respect to the electron coordinates. The derivative is a rank 3 tensor.
|
||||
The first dimension has a dimension of 4 of which the first three coordinates
|
||||
contains the gradient vector and the last index is the laplacian.
|
||||
|
||||
|
||||
\[
|
||||
C_{ij} = \left( 1 - \exp{-\kappa C_{ij}}\right)/\kappa
|
||||
\]
|
||||
|
||||
Here the gradient is defined as follows:
|
||||
|
||||
\[
|
||||
\nabla (C_{ij}(\mathbf{r}_{ee})) = \left(\frac{\delta C_{ij}(\mathbf{r}_{ee})}{\delta x},\frac{\delta C_{ij}(\mathbf{r}_{ee})}{\delta y},\frac{\delta C_{ij}(\mathbf{r}_{ee})}{\delta z} \right)
|
||||
\]
|
||||
and the laplacian is defined as follows:
|
||||
|
||||
\[
|
||||
\triangle (C_{ij}(r_{ee})) = \frac{\delta^2}{\delta x^2} + \frac{\delta^2}{\delta y^2} + \frac{\delta^2}{\delta z^2}
|
||||
\]
|
||||
|
||||
Using the above three formulae, the expression for the gradient and laplacian is
|
||||
as follows:
|
||||
|
||||
\[
|
||||
\frac{\delta C_{ij}(\mathbf{r}_{ee})}{\delta x} = \frac{|(x_i - x_j)|}{r_{ij}} (1 - \kappa R_{ij})
|
||||
\]
|
||||
|
||||
\[
|
||||
\frac{\delta C_{ij}(\mathbf{r}_{ee})}{\delta y} = \frac{|(y_i - y_j)|}{r_{ij}} (1 - \kappa R_{ij})
|
||||
\]
|
||||
|
||||
\[
|
||||
\frac{\delta C_{ij}(\mathbf{r}_{ee})}{\delta z} = \frac{|(z_i - z_j)|}{r_{ij}} (1 - \kappa R_{ij})
|
||||
\]
|
||||
|
||||
\[
|
||||
\Delta(C_{ij}(r_{ee}) = \left[ \frac{2}{r_{ij}} - \kappa \right] (1-\kappa R_{ij})
|
||||
\]
|
||||
|
||||
If the input array is normal (~'N'~), the xyz coordinates are in
|
||||
the leading dimension: ~[n][3]~ in C and ~(3,n)~ in Fortran.
|
||||
|
||||
#+NAME: qmckl_distance_rescaled_deriv_e_args
|
||||
| qmckl_context | context | in | Global state |
|
||||
| char | transa | in | Array ~A~ is ~'N'~: Normal, ~'T'~: Transposed |
|
||||
| char | transb | in | Array ~B~ is ~'N'~: Normal, ~'T'~: Transposed |
|
||||
| int64_t | m | in | Number of points in the first set |
|
||||
| int64_t | n | in | Number of points in the second set |
|
||||
| double | A[][lda] | in | Array containing the $m \times 3$ matrix $A$ |
|
||||
| int64_t | lda | in | Leading dimension of array ~A~ |
|
||||
| double | B[][ldb] | in | Array containing the $n \times 3$ matrix $B$ |
|
||||
| int64_t | ldb | in | Leading dimension of array ~B~ |
|
||||
| double | C[4][n][ldc] | out | Array containing the $4 \times m \times n$ matrix $C$ |
|
||||
| int64_t | ldc | in | Leading dimension of array ~C~ |
|
||||
| double | rescale_factor_kappa | in | Factor for calculating rescaled distances derivatives |
|
||||
|
||||
*** Requirements
|
||||
|
||||
- ~context~ is not ~QMCKL_NULL_CONTEXT~
|
||||
- ~m > 0~
|
||||
- ~n > 0~
|
||||
- ~lda >= 3~ if ~transa == 'N'~
|
||||
- ~lda >= m~ if ~transa == 'T'~
|
||||
- ~ldb >= 3~ if ~transb == 'N'~
|
||||
- ~ldb >= n~ if ~transb == 'T'~
|
||||
- ~ldc >= m~
|
||||
- ~A~ is allocated with at least $3 \times m \times 8$ bytes
|
||||
- ~B~ is allocated with at least $3 \times n \times 8$ bytes
|
||||
- ~C~ is allocated with at least $4 \times m \times n \times 8$ bytes
|
||||
|
||||
*** C header
|
||||
|
||||
#+CALL: generate_c_header(table=qmckl_distance_rescaled_deriv_e_args,rettyp=get_value("CRetType"),fname=get_value("Name"))
|
||||
|
||||
#+RESULTS:
|
||||
#+begin_src c :tangle (eval h_func) :comments org
|
||||
qmckl_exit_code qmckl_distance_rescaled_deriv_e (
|
||||
const qmckl_context context,
|
||||
const char transa,
|
||||
const char transb,
|
||||
const int64_t m,
|
||||
const int64_t n,
|
||||
const double* A,
|
||||
const int64_t lda,
|
||||
const double* B,
|
||||
const int64_t ldb,
|
||||
double* const C,
|
||||
const int64_t ldc,
|
||||
const double rescale_factor_kappa);
|
||||
#+end_src
|
||||
|
||||
*** Source
|
||||
#+begin_src f90 :tangle (eval f)
|
||||
integer function qmckl_distance_rescaled_deriv_e_f(context, transa, transb, m, n, &
|
||||
A, LDA, B, LDB, C, LDC, rescale_factor_kappa) &
|
||||
result(info)
|
||||
use qmckl
|
||||
implicit none
|
||||
integer(qmckl_context) , intent(in) :: context
|
||||
character , intent(in) :: transa, transb
|
||||
integer*8 , intent(in) :: m, n
|
||||
integer*8 , intent(in) :: lda
|
||||
real*8 , intent(in) :: A(lda,*)
|
||||
integer*8 , intent(in) :: ldb
|
||||
real*8 , intent(in) :: B(ldb,*)
|
||||
integer*8 , intent(in) :: ldc
|
||||
real*8 , intent(out) :: C(4,ldc,*)
|
||||
real*8 , intent(in) :: rescale_factor_kappa
|
||||
|
||||
integer*8 :: i,j
|
||||
real*8 :: x, y, z, dist, dist_inv
|
||||
real*8 :: rescale_factor_kappa_inv, rij
|
||||
integer :: transab
|
||||
|
||||
rescale_factor_kappa_inv = 1.0d0/rescale_factor_kappa;
|
||||
|
||||
info = QMCKL_SUCCESS
|
||||
|
||||
if (context == QMCKL_NULL_CONTEXT) then
|
||||
info = QMCKL_INVALID_CONTEXT
|
||||
return
|
||||
endif
|
||||
|
||||
if (m <= 0_8) then
|
||||
info = QMCKL_INVALID_ARG_4
|
||||
return
|
||||
endif
|
||||
|
||||
if (n <= 0_8) then
|
||||
info = QMCKL_INVALID_ARG_5
|
||||
return
|
||||
endif
|
||||
|
||||
if (transa == 'N' .or. transa == 'n') then
|
||||
transab = 0
|
||||
else if (transa == 'T' .or. transa == 't') then
|
||||
transab = 1
|
||||
else
|
||||
transab = -100
|
||||
endif
|
||||
|
||||
if (transb == 'N' .or. transb == 'n') then
|
||||
continue
|
||||
else if (transa == 'T' .or. transa == 't') then
|
||||
transab = transab + 2
|
||||
else
|
||||
transab = -100
|
||||
endif
|
||||
|
||||
! check for LDA
|
||||
if (transab < 0) then
|
||||
info = QMCKL_INVALID_ARG_1
|
||||
return
|
||||
endif
|
||||
|
||||
if (iand(transab,1) == 0 .and. LDA < 3) then
|
||||
info = QMCKL_INVALID_ARG_7
|
||||
return
|
||||
endif
|
||||
|
||||
if (iand(transab,1) == 1 .and. LDA < m) then
|
||||
info = QMCKL_INVALID_ARG_7
|
||||
return
|
||||
endif
|
||||
|
||||
if (iand(transab,2) == 0 .and. LDA < 3) then
|
||||
info = QMCKL_INVALID_ARG_7
|
||||
return
|
||||
endif
|
||||
|
||||
if (iand(transab,2) == 2 .and. LDA < m) then
|
||||
info = QMCKL_INVALID_ARG_7
|
||||
return
|
||||
endif
|
||||
|
||||
! check for LDB
|
||||
if (iand(transab,1) == 0 .and. LDB < 3) then
|
||||
info = QMCKL_INVALID_ARG_9
|
||||
return
|
||||
endif
|
||||
|
||||
if (iand(transab,1) == 1 .and. LDB < n) then
|
||||
info = QMCKL_INVALID_ARG_9
|
||||
return
|
||||
endif
|
||||
|
||||
if (iand(transab,2) == 0 .and. LDB < 3) then
|
||||
info = QMCKL_INVALID_ARG_9
|
||||
return
|
||||
endif
|
||||
|
||||
if (iand(transab,2) == 2 .and. LDB < n) then
|
||||
info = QMCKL_INVALID_ARG_9
|
||||
return
|
||||
endif
|
||||
|
||||
! check for LDC
|
||||
if (LDC < m) then
|
||||
info = QMCKL_INVALID_ARG_11
|
||||
return
|
||||
endif
|
||||
|
||||
|
||||
select case (transab)
|
||||
|
||||
case(0)
|
||||
|
||||
do j=1,n
|
||||
do i=1,m
|
||||
x = A(1,i) - B(1,j)
|
||||
y = A(2,i) - B(2,j)
|
||||
z = A(3,i) - B(3,j)
|
||||
dist = dsqrt(x*x + y*y + z*z)
|
||||
dist_inv = 1.0d0/dist
|
||||
rij = (1.0d0 - dexp(-rescale_factor_kappa * dist)) * rescale_factor_kappa_inv
|
||||
C(1,i,j) = x * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(2,i,j) = y * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(3,i,j) = z * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(4,i,j) = (2.0d0 * dist_inv - rescale_factor_kappa_inv) * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
end do
|
||||
end do
|
||||
|
||||
case(1)
|
||||
|
||||
do j=1,n
|
||||
do i=1,m
|
||||
x = A(i,1) - B(1,j)
|
||||
y = A(i,2) - B(2,j)
|
||||
z = A(i,3) - B(3,j)
|
||||
dist = dsqrt(x*x + y*y + z*z)
|
||||
dist_inv = 1.0d0/dist
|
||||
rij = (1.0d0 - dexp(-rescale_factor_kappa * dist)) * rescale_factor_kappa_inv
|
||||
C(1,i,j) = x * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(2,i,j) = y * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(3,i,j) = z * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(4,i,j) = (2.0d0 * dist_inv - rescale_factor_kappa_inv) * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
end do
|
||||
end do
|
||||
|
||||
case(2)
|
||||
|
||||
do j=1,n
|
||||
do i=1,m
|
||||
x = A(1,i) - B(j,1)
|
||||
y = A(2,i) - B(j,2)
|
||||
z = A(3,i) - B(j,3)
|
||||
dist = dsqrt(x*x + y*y + z*z)
|
||||
dist_inv = 1.0d0/dist
|
||||
rij = (1.0d0 - dexp(-rescale_factor_kappa * dist)) * rescale_factor_kappa_inv
|
||||
C(1,i,j) = x * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(2,i,j) = y * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(3,i,j) = z * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(4,i,j) = (2.0d0 * dist_inv - rescale_factor_kappa_inv) * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
end do
|
||||
end do
|
||||
|
||||
case(3)
|
||||
|
||||
do j=1,n
|
||||
do i=1,m
|
||||
x = A(i,1) - B(j,1)
|
||||
y = A(i,2) - B(j,2)
|
||||
z = A(i,3) - B(j,3)
|
||||
dist = dsqrt(x*x + y*y + z*z)
|
||||
dist_inv = 1.0d0/dist
|
||||
rij = (1.0d0 - dexp(-rescale_factor_kappa * dist)) * rescale_factor_kappa_inv
|
||||
C(1,i,j) = x * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(2,i,j) = y * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(3,i,j) = z * dist_inv * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
C(4,i,j) = (2.0d0 * dist_inv - rescale_factor_kappa_inv) * ( 1.0d0 - rescale_factor_kappa_inv * rij)
|
||||
end do
|
||||
end do
|
||||
|
||||
end select
|
||||
|
||||
end function qmckl_distance_rescaled_deriv_e_f
|
||||
#+end_src
|
||||
|
||||
*** Performance
|
||||
|
||||
This function is more efficient when ~A~ and ~B~ are transposed.
|
||||
|
||||
** C interface :noexport:
|
||||
|
||||
#+CALL: generate_c_interface(table=qmckl_distance_rescaled_deriv_e_args,rettyp=get_value("FRetType"),fname=get_value("Name"))
|
||||
|
||||
#+RESULTS:
|
||||
#+begin_src f90 :tangle (eval f) :comments org :exports none
|
||||
integer(c_int32_t) function qmckl_distance_rescaled_deriv_e &
|
||||
(context, transa, transb, m, n, A, lda, B, ldb, C, ldc, rescale_factor_kappa) &
|
||||
bind(C) result(info)
|
||||
|
||||
use, intrinsic :: iso_c_binding
|
||||
implicit none
|
||||
|
||||
integer (c_int64_t) , intent(in) , value :: context
|
||||
character , intent(in) , value :: transa
|
||||
character , intent(in) , value :: transb
|
||||
integer (c_int64_t) , intent(in) , value :: m
|
||||
integer (c_int64_t) , intent(in) , value :: n
|
||||
real (c_double ) , intent(in) :: A(lda,*)
|
||||
integer (c_int64_t) , intent(in) , value :: lda
|
||||
real (c_double ) , intent(in) :: B(ldb,*)
|
||||
integer (c_int64_t) , intent(in) , value :: ldb
|
||||
real (c_double ) , intent(out) :: C(4,ldc,n)
|
||||
integer (c_int64_t) , intent(in) , value :: ldc
|
||||
real (c_double ) , intent(in) , value :: rescale_factor_kappa
|
||||
|
||||
integer(c_int32_t), external :: qmckl_distance_rescaled_deriv_e_f
|
||||
info = qmckl_distance_rescaled_deriv_e_f &
|
||||
(context, transa, transb, m, n, A, lda, B, ldb, C, ldc, rescale_factor_kappa)
|
||||
|
||||
end function qmckl_distance_rescaled_deriv_e
|
||||
#+end_src
|
||||
|
||||
#+CALL: generate_f_interface(table=qmckl_distance_rescaled_deriv_e_args,rettyp=get_value("FRetType"),fname=get_value("Name"))
|
||||
|
||||
#+RESULTS:
|
||||
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
||||
interface
|
||||
integer(c_int32_t) function qmckl_distance_rescaled_deriv_e &
|
||||
(context, transa, transb, m, n, A, lda, B, ldb, C, ldc, rescale_factor_kappa) &
|
||||
bind(C)
|
||||
use, intrinsic :: iso_c_binding
|
||||
import
|
||||
implicit none
|
||||
|
||||
integer (c_int64_t) , intent(in) , value :: context
|
||||
character , intent(in) , value :: transa
|
||||
character , intent(in) , value :: transb
|
||||
integer (c_int64_t) , intent(in) , value :: m
|
||||
integer (c_int64_t) , intent(in) , value :: n
|
||||
real (c_double ) , intent(in) :: A(lda,*)
|
||||
integer (c_int64_t) , intent(in) , value :: lda
|
||||
real (c_double ) , intent(in) :: B(ldb,*)
|
||||
integer (c_int64_t) , intent(in) , value :: ldb
|
||||
real (c_double ) , intent(out) :: C(4,ldc,n)
|
||||
integer (c_int64_t) , intent(in) , value :: ldc
|
||||
real (c_double ) , intent(in) , value :: rescale_factor_kappa
|
||||
|
||||
end function qmckl_distance_rescaled_deriv_e
|
||||
end interface
|
||||
#+end_src
|
||||
|
||||
* End of files :noexport:
|
||||
|
||||
#+begin_src c :comments link :tangle (eval c_test)
|
||||
|
Loading…
Reference in New Issue
Block a user