-
Notifications
You must be signed in to change notification settings - Fork 2
/
el_nuc_el_blas.irp.f
104 lines (86 loc) · 3.19 KB
/
el_nuc_el_blas.irp.f
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
BEGIN_PROVIDER [ double precision, tmp_c, (nelec_8,nnuc,0:ncord,0:ncord-1) ]
&BEGIN_PROVIDER [ double precision, dtmp_c, (nelec_8,4,nnuc,0:ncord,0:ncord-1) ]
implicit none
BEGIN_DOC
! Calculate the intermediate buffers
! tmp_c:
! r_{ij}^k . R_{ja}^l -> tmp_c_{ia}^{kl}
!
! dtmp_c:
! dr_{ij}^k . R_{ja}^l -> dtmp_c_{ia}^{kl}
END_DOC
integer :: k
! r_{ij}^k . R_{ja}^l -> tmp_c_{ia}^{kl}
do k=0,ncord-1
call dgemm('N','N', nelec, nnuc*(ncord+1), nelec, 1.d0, &
rescale_een_e(1,1,k), size(rescale_een_e,1), &
rescale_een_n(1,1,0), size(rescale_een_n,1), 0.d0, &
tmp_c(1,1,0,k), size(tmp_c,1))
enddo
! dr_{ij}^k . R_{ja}^l -> dtmp_c_{ia}^{kl}
do k=0,ncord-1
call dgemm('N','N', 4*nelec_8, nnuc*(ncord+1), nelec, 1.d0, &
rescale_een_e_deriv_e(1,1,1,k), &
size(rescale_een_e_deriv_e,1)*size(rescale_een_e_deriv_e,2), &
rescale_een_n(1,1,0), &
size(rescale_een_n,1), 0.d0, &
dtmp_c(1,1,1,0,k), size(dtmp_c,1)*size(dtmp_c,2))
enddo
END_PROVIDER
BEGIN_PROVIDER [ double precision, factor_een_blas ]
&BEGIN_PROVIDER [ double precision, factor_een_deriv_e_blas, (nelec_8,4) ]
implicit none
BEGIN_DOC
! Dimensions 1-3 : dx, dy, dz
! Dimension 4 : d2x + d2y + d2z
END_DOC
integer :: i, j, a, p, k, l, lmax, m, n, ii
double precision :: accu, cn, cn2
! double precision,dimension(:),allocatable :: cn
factor_een_blas = 0.0d0
factor_een_deriv_e_blas(:,:) = 0.0d0
do n = 1, dim_cord_vect
l = lkpm_of_cindex(1,n)
k = lkpm_of_cindex(2,n)
p = lkpm_of_cindex(3,n)
m = lkpm_of_cindex(4,n)
do a = 1, nnuc
cn = cord_vect_full(n, a)
if (cn == 0.d0) cycle
accu = 0.d0
do j=1,nelec
accu = accu + rescale_een_n(j,a,m) * tmp_c(j,a,m+l,k)
enddo
factor_een_blas = factor_een_blas + accu * cn
cn2 = cn+cn
do ii=1,4
do j=1,nelec
factor_een_deriv_e_blas(j,ii) = factor_een_deriv_e_blas(j,ii) + (&
tmp_c (j,a,m+l,k) *rescale_een_n_deriv_e(j,ii,a,m) + &
dtmp_c(j,ii,a,m+l,k) * rescale_een_n(j,a,m) ) * cn
enddo
enddo
do j=1,nelec
factor_een_deriv_e_blas(j,4) = factor_een_deriv_e_blas(j,4) + (&
dtmp_c(j,1,a,m+l,k) * rescale_een_n_deriv_e(j,1,a,m ) + &
dtmp_c(j,2,a,m+l,k) * rescale_een_n_deriv_e(j,2,a,m ) + &
dtmp_c(j,3,a,m+l,k) * rescale_een_n_deriv_e(j,3,a,m ) &
)*cn2
enddo
do ii=1,4
do j=1,nelec
factor_een_deriv_e_blas(j,ii) = factor_een_deriv_e_blas(j,ii) + (&
tmp_c(j,a,m,k) * rescale_een_n_deriv_e(j,ii,a,m+l) + &
dtmp_c(j,ii,a,m,k) * rescale_een_n(j,a,m+l) ) * cn
enddo
enddo
do j=1,nelec
factor_een_deriv_e_blas(j,4) = factor_een_deriv_e_blas(j,4) + (&
dtmp_c(j,1,a,m ,k) * rescale_een_n_deriv_e(j,1,a,m+l) + &
dtmp_c(j,2,a,m ,k) * rescale_een_n_deriv_e(j,2,a,m+l) + &
dtmp_c(j,3,a,m ,k) * rescale_een_n_deriv_e(j,3,a,m+l) &
)*cn2
enddo
enddo
enddo
END_PROVIDER