(git:374b731)
Loading...
Searching...
No Matches
lri_integrals.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Calculates integrals for LRIGPW method
10!> lri : local resolution of the identity
11!> \par History
12!> created JGH [08.2012]
13!> Dorothea Golze [02.2014] (1) extended, re-structured, cleaned
14!> (2) heavily debugged
15!> split off JGH [11.2017]
16!> \authors JGH
17!> Dorothea Golze
18! **************************************************************************************************
19MODULE lri_integrals
20 USE basis_set_types, ONLY: gto_basis_set_type
21 USE generic_os_integrals, ONLY: int_overlap_ab_os,&
22 int_overlap_aba_os,&
23 int_overlap_abb_os
24 USE generic_shg_integrals, ONLY: get_abb_same_kind,&
25 int_overlap_ab_shg_low,&
26 int_overlap_aba_shg_low,&
27 int_overlap_abb_shg_low,&
28 lri_precalc_angular_shg_part
29 USE kinds, ONLY: dp
30 USE lri_environment_types, ONLY: lri_environment_type,&
31 lri_int_type
32#include "./base/base_uses.f90"
33
34 IMPLICIT NONE
35
36 PRIVATE
37
38! **************************************************************************************************
39 TYPE int_type
40 REAL(kind=dp), DIMENSION(:, :), ALLOCATABLE :: sabint
41 REAL(kind=dp), DIMENSION(:, :), ALLOCATABLE :: sooint
42 REAL(kind=dp), DIMENSION(:, :, :), ALLOCATABLE :: abaint
43 REAL(kind=dp), DIMENSION(:, :, :), ALLOCATABLE :: abbint
44 END TYPE int_type
45 !
46 TYPE dint_type
47 REAL(kind=dp), DIMENSION(:, :, :), ALLOCATABLE :: dsabint
48 REAL(kind=dp), DIMENSION(:, :, :), ALLOCATABLE :: dsooint
49 REAL(kind=dp), DIMENSION(:, :, :, :), ALLOCATABLE :: dabaint
50 REAL(kind=dp), DIMENSION(:, :, :, :), ALLOCATABLE :: dabbint
51 END TYPE dint_type
52! **************************************************************************************************
53
54 CHARACTER(len=*), PARAMETER, PRIVATE :: modulen = 'lri_integrals'
55
56 PUBLIC :: lri_int, lri_int2, lri_dint, lri_dint2, int_type, dint_type, allocate_int_type, deallocate_int_type
57
58! **************************************************************************************************
59
60CONTAINS
61
62! **************************************************************************************************
63!> \brief calcuates the lri integrals using solid harmonic Gaussians
64!> \param lri_env ...
65!> \param lrii ...
66!> \param rab distance vector
67!> \param obasa orb basis on A
68!> \param obasb orb basis on B
69!> \param fbasa aux basis on A
70!> \param fbasb aux basis on B
71!> \param iatom index atom A
72!> \param jatom index atom B
73!> \param ikind kind atom A
74!> \param jkind kind atom B
75!> \param calculate_forces ...
76! **************************************************************************************************
77 SUBROUTINE lri_int(lri_env, lrii, rab, obasa, obasb, fbasa, fbasb, &
78 iatom, jatom, ikind, jkind, calculate_forces)
79
80 TYPE(lri_environment_type), POINTER :: lri_env
81 TYPE(lri_int_type), POINTER :: lrii
82 REAL(kind=dp), DIMENSION(3), INTENT(IN) :: rab
83 TYPE(gto_basis_set_type), POINTER :: obasa, obasb, fbasa, fbasb
84 INTEGER, INTENT(IN) :: iatom, jatom, ikind, jkind
85 LOGICAL, INTENT(IN) :: calculate_forces
86
87 IF (lri_env%use_shg_integrals) THEN
88 CALL lri_int_shg(lri_env, lrii, rab, obasa, obasb, fbasa, fbasb, &
89 iatom, jatom, ikind, jkind, calculate_forces)
90 ELSE
91 CALL lri_int_os(lri_env, lrii, rab, obasa, obasb, fbasa, fbasb, &
92 iatom, jatom, ikind, calculate_forces)
93 END IF
94
95 END SUBROUTINE lri_int
96
97! **************************************************************************************************
98!> \brief calcuates the lri integrals using solid harmonic Gaussians
99!> \param lri_env ...
100!> \param lrii ...
101!> \param lriint ...
102!> \param rab distance vector
103!> \param obasa orb basis on A
104!> \param obasb orb basis on B
105!> \param fbasa aux basis on A
106!> \param fbasb aux basis on B
107!> \param iatom index atom A
108!> \param jatom index atom B
109!> \param ikind kind atom A
110!> \param jkind kind atom B
111! **************************************************************************************************
112 SUBROUTINE lri_int2(lri_env, lrii, lriint, rab, obasa, obasb, fbasa, fbasb, iatom, jatom, ikind, jkind)
113
114 TYPE(lri_environment_type), POINTER :: lri_env
115 TYPE(lri_int_type), POINTER :: lrii
116 TYPE(int_type) :: lriint
117 REAL(kind=dp), DIMENSION(3), INTENT(IN) :: rab
118 TYPE(gto_basis_set_type), POINTER :: obasa, obasb, fbasa, fbasb
119 INTEGER, INTENT(IN) :: iatom, jatom, ikind, jkind
120
121 INTEGER :: nba, nfa
122 INTEGER, DIMENSION(:, :, :), POINTER :: fba_index, fbb_index, oba_index, &
123 obb_index
124 REAL(kind=dp) :: dab
125 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :) :: dummy1, waux_mat
126 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :, :) :: dummy2, dwaux_mat
127 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: scon_fba, scon_fbb, scon_oba, scon_obb
128 REAL(kind=dp), DIMENSION(:, :, :, :), POINTER :: scona_mix, sconb_mix
129
130 dab = sqrt(sum(rab*rab))
131 IF (iatom == jatom .AND. dab < lri_env%delta) THEN
132 cpassert(ALLOCATED(lriint%sabint))
133 cpassert(ALLOCATED(lriint%sooint))
134 cpassert(ALLOCATED(lriint%abaint))
135 nba = obasa%nsgf
136 nfa = fbasa%nsgf
137 lriint%sabint(1:nfa, 1:nfa) = lri_env%bas_prop(ikind)%ri_ovlp(1:nfa, 1:nfa)
138 lriint%sooint(1:nba, 1:nba) = lri_env%bas_prop(ikind)%orb_ovlp(1:nba, 1:nba)
139 lriint%abaint(1:nba, 1:nba, 1:nfa) = lri_env%bas_prop(ikind)%ovlp3
140 ELSE
141 IF (lri_env%use_shg_integrals) THEN
142 scon_oba => lri_env%bas_prop(ikind)%scon_orb
143 scon_obb => lri_env%bas_prop(jkind)%scon_orb
144 scon_fba => lri_env%bas_prop(ikind)%scon_ri
145 scon_fbb => lri_env%bas_prop(jkind)%scon_ri
146 scona_mix => lri_env%bas_prop(ikind)%scon_mix
147 sconb_mix => lri_env%bas_prop(jkind)%scon_mix
148 oba_index => lri_env%bas_prop(ikind)%orb_index
149 fba_index => lri_env%bas_prop(ikind)%ri_index
150 obb_index => lri_env%bas_prop(jkind)%orb_index
151 fbb_index => lri_env%bas_prop(jkind)%ri_index
152 !
153 CALL lri_precalc_angular_shg_part(obasa, obasb, fbasa, fbasb, rab, waux_mat, dwaux_mat, &
154 .false.)
155 !*** (fa,fb)
156 IF (ALLOCATED(lriint%sabint)) THEN
157 CALL int_overlap_ab_shg_low(lriint%sabint, dummy1, rab, fbasa, fbasb, scon_fba, scon_fbb, &
158 waux_mat, dwaux_mat, .true., .false., contraction_high=.false.)
159 END IF
160 !*** (a,b)
161 IF (ALLOCATED(lriint%sooint)) THEN
162 CALL int_overlap_ab_shg_low(lriint%sooint, dummy1, rab, obasa, obasb, scon_oba, scon_obb, &
163 waux_mat, dwaux_mat, .true., .false., contraction_high=.true.)
164 END IF
165 !*** (a,b,fa)
166 IF (ALLOCATED(lriint%abaint)) THEN
167 CALL int_overlap_aba_shg_low(lriint%abaint, dummy2, rab, obasa, obasb, fbasa, &
168 scon_obb, scona_mix, oba_index, fba_index, &
169 lri_env%cg_shg%cg_coeff, lri_env%cg_shg%cg_none0_list, &
170 lri_env%cg_shg%ncg_none0, &
171 waux_mat, dwaux_mat, .true., .false.)
172 END IF
173 !*** (a,b,fb)
174 IF (ALLOCATED(lriint%abbint)) THEN
175 IF (ikind == jkind) THEN
176 cpassert(ALLOCATED(lriint%abaint))
177 CALL get_abb_same_kind(lriint%abbint, dummy2, lriint%abaint, dummy2, &
178 rab, obasa, fbasa, .true., .false.)
179 ELSE
180 CALL int_overlap_abb_shg_low(lriint%abbint, dummy2, rab, obasa, obasb, fbasb, &
181 scon_oba, sconb_mix, obb_index, fbb_index, &
182 lri_env%cg_shg%cg_coeff, lri_env%cg_shg%cg_none0_list, &
183 lri_env%cg_shg%ncg_none0, &
184 waux_mat, dwaux_mat, .true., .false.)
185 END IF
186 END IF
187 DEALLOCATE (waux_mat, dwaux_mat)
188 ELSE
189 !*** (fa,fb)
190 IF (ALLOCATED(lriint%sabint)) THEN
191 CALL int_overlap_ab_os(lriint%sabint, dummy1, rab, fbasa, fbasb, &
192 .false., lri_env%debug, lrii%dmax_ab)
193 END IF
194 !*** (a,b)
195 IF (ALLOCATED(lriint%sooint)) THEN
196 CALL int_overlap_ab_os(lriint%sooint, dummy1, rab, obasa, obasb, &
197 .false., lri_env%debug, lrii%dmax_oo)
198 END IF
199 !*** (a,b,fa)
200 IF (ALLOCATED(lriint%abaint)) THEN
201 CALL int_overlap_aba_os(lriint%abaint, dummy2, rab, obasa, obasb, fbasa, &
202 .false., lri_env%debug, lrii%dmax_aba)
203 END IF
204 !*** (a,b,fb)
205 IF (ALLOCATED(lriint%abbint)) THEN
206 CALL int_overlap_abb_os(lriint%abbint, dummy2, rab, obasa, obasb, fbasb, &
207 .false., lri_env%debug, lrii%dmax_abb)
208 END IF
209 END IF
210 END IF
211
212 END SUBROUTINE lri_int2
213
214! **************************************************************************************************
215
216! **************************************************************************************************
217!> \brief ...
218!> \param lri_env ...
219!> \param lrii ...
220!> \param rab ...
221!> \param obasa ...
222!> \param obasb ...
223!> \param fbasa ...
224!> \param fbasb ...
225!> \param iatom ...
226!> \param jatom ...
227!> \param ikind ...
228!> \param jkind ...
229!> \param calculate_forces ...
230! **************************************************************************************************
231 SUBROUTINE lri_dint(lri_env, lrii, rab, obasa, obasb, fbasa, fbasb, &
232 iatom, jatom, ikind, jkind, calculate_forces)
233
234 TYPE(lri_environment_type), POINTER :: lri_env
235 TYPE(lri_int_type), POINTER :: lrii
236 REAL(kind=dp), DIMENSION(3), INTENT(IN) :: rab
237 TYPE(gto_basis_set_type), POINTER :: obasa, obasb, fbasa, fbasb
238 INTEGER, INTENT(IN) :: iatom, jatom, ikind, jkind
239 LOGICAL, INTENT(IN) :: calculate_forces
240
241 IF (lri_env%use_shg_integrals) THEN
242 CALL lri_int_shg(lri_env, lrii, rab, obasa, obasb, fbasa, fbasb, &
243 iatom, jatom, ikind, jkind, calculate_forces)
244 ELSE
245 CALL lri_int_os(lri_env, lrii, rab, obasa, obasb, fbasa, fbasb, &
246 iatom, jatom, ikind, calculate_forces)
247 END IF
248
249 END SUBROUTINE lri_dint
250
251! **************************************************************************************************
252!> \brief ...
253!> \param lri_env ...
254!> \param lrii ...
255!> \param lridint ...
256!> \param rab ...
257!> \param obasa ...
258!> \param obasb ...
259!> \param fbasa ...
260!> \param fbasb ...
261!> \param iatom ...
262!> \param jatom ...
263!> \param ikind ...
264!> \param jkind ...
265! **************************************************************************************************
266 SUBROUTINE lri_dint2(lri_env, lrii, lridint, rab, obasa, obasb, fbasa, fbasb, &
267 iatom, jatom, ikind, jkind)
268
269 TYPE(lri_environment_type), POINTER :: lri_env
270 TYPE(lri_int_type), POINTER :: lrii
271 TYPE(dint_type) :: lridint
272 REAL(kind=dp), DIMENSION(3), INTENT(IN) :: rab
273 TYPE(gto_basis_set_type), POINTER :: obasa, obasb, fbasa, fbasb
274 INTEGER, INTENT(IN) :: iatom, jatom, ikind, jkind
275
276 INTEGER :: nba, nbb, nfa, nfb
277 INTEGER, DIMENSION(:, :, :), POINTER :: fba_index, fbb_index, oba_index, &
278 obb_index
279 REAL(kind=dp) :: dab
280 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: dummy1
281 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :) :: dummy2, waux_mat
282 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :, :) :: dwaux_mat
283 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: scon_fba, scon_fbb, scon_oba, scon_obb
284 REAL(kind=dp), DIMENSION(:, :, :, :), POINTER :: scona_mix, sconb_mix
285
286 dab = sqrt(sum(rab*rab))
287 nba = obasa%nsgf
288 nbb = obasb%nsgf
289 nfa = fbasa%nsgf
290 nfb = fbasb%nsgf
291
292 IF (iatom == jatom .AND. dab < lri_env%delta) THEN
293 ! nothing to do
294 ELSE
295 IF (lrii%calc_force_pair) THEN
296 IF (lri_env%use_shg_integrals) THEN
297 scon_oba => lri_env%bas_prop(ikind)%scon_orb
298 scon_obb => lri_env%bas_prop(jkind)%scon_orb
299 scon_fba => lri_env%bas_prop(ikind)%scon_ri
300 scon_fbb => lri_env%bas_prop(jkind)%scon_ri
301 scona_mix => lri_env%bas_prop(ikind)%scon_mix
302 sconb_mix => lri_env%bas_prop(jkind)%scon_mix
303 oba_index => lri_env%bas_prop(ikind)%orb_index
304 fba_index => lri_env%bas_prop(ikind)%ri_index
305 obb_index => lri_env%bas_prop(jkind)%orb_index
306 fbb_index => lri_env%bas_prop(jkind)%ri_index
307
308 CALL lri_precalc_angular_shg_part(obasa, obasb, fbasa, fbasb, rab, &
309 waux_mat, dwaux_mat, .true.)
310 !*** (fa,fb)
311 IF (ALLOCATED(lridint%dsabint)) THEN
312 CALL int_overlap_ab_shg_low(dummy1, lridint%dsabint, rab, fbasa, fbasb, scon_fba, scon_fbb, &
313 waux_mat, dwaux_mat, .false., .true., contraction_high=.false.)
314 END IF
315 !*** (a,b)
316 IF (ALLOCATED(lridint%dsooint)) THEN
317 CALL int_overlap_ab_shg_low(dummy1, lridint%dsooint, rab, obasa, obasb, scon_oba, scon_obb, &
318 waux_mat, dwaux_mat, .false., .true., contraction_high=.true.)
319 END IF
320 !*** (a,b,fa)
321 IF (ALLOCATED(lridint%dabaint)) THEN
322 CALL int_overlap_aba_shg_low(dummy2, lridint%dabaint, rab, obasa, obasb, fbasa, &
323 scon_obb, scona_mix, oba_index, fba_index, &
324 lri_env%cg_shg%cg_coeff, lri_env%cg_shg%cg_none0_list, &
325 lri_env%cg_shg%ncg_none0, &
326 waux_mat, dwaux_mat, .false., .true.)
327 END IF
328 !*** (a,b,fb)
329 IF (ALLOCATED(lridint%dabbint)) THEN
330 IF (ikind == jkind) THEN
331 cpassert(ALLOCATED(lridint%dabaint))
332 CALL get_abb_same_kind(dummy2, lridint%dabbint, dummy2, lridint%dabaint, &
333 rab, obasa, fbasa, .false., .true.)
334 ELSE
335 CALL int_overlap_abb_shg_low(dummy2, lridint%dabbint, rab, obasa, obasb, fbasb, &
336 scon_oba, sconb_mix, obb_index, fbb_index, &
337 lri_env%cg_shg%cg_coeff, lri_env%cg_shg%cg_none0_list, &
338 lri_env%cg_shg%ncg_none0, &
339 waux_mat, dwaux_mat, .false., .true.)
340 END IF
341 END IF
342 DEALLOCATE (waux_mat, dwaux_mat)
343
344 ELSE
345
346 !*** (fa,fb)
347 IF (ALLOCATED(lridint%dsabint)) THEN
348 ALLOCATE (dummy1(nfa, nfb))
349 CALL int_overlap_ab_os(dummy1, lridint%dsabint, rab, fbasa, fbasb, &
350 .true., lri_env%debug, lrii%dmax_ab)
351 DEALLOCATE (dummy1)
352 END IF
353 !*** (a,b)
354 IF (ALLOCATED(lridint%dsooint)) THEN
355 ALLOCATE (dummy1(nba, nbb))
356 CALL int_overlap_ab_os(dummy1, lridint%dsooint, rab, obasa, obasb, &
357 .true., lri_env%debug, lrii%dmax_oo)
358 DEALLOCATE (dummy1)
359 END IF
360 !*** (a,b,fa)
361 IF (ALLOCATED(lridint%dabaint)) THEN
362 ALLOCATE (dummy2(nba, nbb, nfa))
363 CALL int_overlap_aba_os(dummy2, lridint%dabaint, rab, obasa, obasb, fbasa, &
364 .true., lri_env%debug, lrii%dmax_aba)
365 DEALLOCATE (dummy2)
366 END IF
367 !*** (a,b,fb)
368 IF (ALLOCATED(lridint%dabbint)) THEN
369 ALLOCATE (dummy2(nba, nbb, nfb))
370 CALL int_overlap_abb_os(dummy2, lridint%dabbint, rab, obasa, obasb, fbasb, &
371 .true., lri_env%debug, lrii%dmax_abb)
372 DEALLOCATE (dummy2)
373 END IF
374 END IF
375 END IF
376 END IF
377
378 END SUBROUTINE lri_dint2
379
380! **************************************************************************************************
381!> \brief calculates the lri intergrals according to the Obara-Saika (OS)
382!> scheme
383!> \param lri_env ...
384!> \param lrii ...
385!> \param rab distance vector
386!> \param obasa orb basis on center A
387!> \param obasb orb basis on center B
388!> \param fbasa aux basis on center A
389!> \param fbasb aux basis on center B
390!> \param iatom index atom A
391!> \param jatom index atom B
392!> \param ikind kind atom A
393!> \param calculate_forces ...
394! **************************************************************************************************
395 SUBROUTINE lri_int_os(lri_env, lrii, rab, obasa, obasb, fbasa, fbasb, &
396 iatom, jatom, ikind, calculate_forces)
397
398 TYPE(lri_environment_type), POINTER :: lri_env
399 TYPE(lri_int_type), POINTER :: lrii
400 REAL(kind=dp), DIMENSION(3), INTENT(IN) :: rab
401 TYPE(gto_basis_set_type), POINTER :: obasa, obasb, fbasa, fbasb
402 INTEGER, INTENT(IN) :: iatom, jatom, ikind
403 LOGICAL, INTENT(IN) :: calculate_forces
404
405 CHARACTER(LEN=*), PARAMETER :: routinen = 'lri_int_os'
406
407 INTEGER :: handle, nba, nbb, nfa, nfb
408 LOGICAL :: do_force
409 REAL(kind=dp) :: dab
410
411 CALL timeset(routinen, handle)
412
413 dab = sqrt(sum(rab*rab))
414 nba = obasa%nsgf
415 nbb = obasb%nsgf
416 nfa = fbasa%nsgf
417 nfb = fbasb%nsgf
418
419 !*** calculate overlap integrals; for iatom=jatom this is the self-overlap
420 IF (iatom == jatom .AND. dab < lri_env%delta) THEN
421 !*** (fa,fa)
422 lrii%sab(1:nfa, 1:nfa) = lri_env%bas_prop(ikind)%ri_ovlp(1:nfa, 1:nfa)
423 !*** (a,a)
424 lrii%soo(1:nba, 1:nba) = lri_env%bas_prop(ikind)%orb_ovlp(1:nba, 1:nba)
425 !*** (a,a,fa)
426 CALL int_overlap_aba_os(lrii%abaint, rab=rab, oba=obasa, obb=obasb, &
427 fba=fbasa, calculate_forces=.false., debug=lri_env%debug, &
428 dmax=lrii%dmax_aba)
429 IF (calculate_forces) THEN
430 lrii%dsab = 0._dp
431 lrii%dsoo = 0._dp
432 lrii%dabdaint = 0.0_dp
433 lrii%dabbint = 0.0_dp
434 END IF
435 ELSE
436 IF (calculate_forces) THEN
437 do_force = lrii%calc_force_pair
438 ELSE
439 do_force = .false.
440 END IF
441 !*** (fa,fb)
442 CALL int_overlap_ab_os(lrii%sab, lrii%dsab, rab, fbasa, fbasb, &
443 do_force, lri_env%debug, lrii%dmax_ab)
444 !*** (a,b)
445 CALL int_overlap_ab_os(lrii%soo, lrii%dsoo, rab, obasa, obasb, &
446 do_force, lri_env%debug, lrii%dmax_oo)
447 !*** (a,b,fa)
448 CALL int_overlap_aba_os(lrii%abaint, lrii%dabdaint, rab, obasa, obasb, fbasa, &
449 do_force, lri_env%debug, lrii%dmax_aba)
450 !*** (a,b,fb)
451 CALL int_overlap_abb_os(lrii%abbint, lrii%dabbint, rab, obasa, obasb, fbasb, &
452 do_force, lri_env%debug, lrii%dmax_abb)
453 END IF
454
455 CALL timestop(handle)
456
457 END SUBROUTINE lri_int_os
458
459! **************************************************************************************************
460!> \brief calcuates the lri integrals using solid harmonic Gaussians
461!> \param lri_env ...
462!> \param lrii ...
463!> \param rab distance vector
464!> \param obasa orb basis on A
465!> \param obasb orb basis on B
466!> \param fbasa aux basis on A
467!> \param fbasb aux basis on B
468!> \param iatom index atom A
469!> \param jatom index atom B
470!> \param ikind kind atom A
471!> \param jkind kind atom B
472!> \param calculate_forces ...
473! **************************************************************************************************
474 SUBROUTINE lri_int_shg(lri_env, lrii, rab, obasa, obasb, fbasa, fbasb, &
475 iatom, jatom, ikind, jkind, calculate_forces)
476
477 TYPE(lri_environment_type), POINTER :: lri_env
478 TYPE(lri_int_type), POINTER :: lrii
479 REAL(kind=dp), DIMENSION(3), INTENT(IN) :: rab
480 TYPE(gto_basis_set_type), POINTER :: obasa, obasb, fbasa, fbasb
481 INTEGER, INTENT(IN) :: iatom, jatom, ikind, jkind
482 LOGICAL, INTENT(IN) :: calculate_forces
483
484 CHARACTER(LEN=*), PARAMETER :: routinen = 'lri_int_shg'
485
486 INTEGER :: handle, nba, nbb, nfa, nfb
487 INTEGER, DIMENSION(:, :, :), POINTER :: fba_index, fbb_index, oba_index, &
488 obb_index
489 LOGICAL :: do_force
490 REAL(kind=dp) :: dab
491 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :) :: waux_mat
492 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :, :) :: dwaux_mat
493 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: scon_fba, scon_fbb, scon_oba, scon_obb
494 REAL(kind=dp), DIMENSION(:, :, :, :), POINTER :: scona_mix, sconb_mix
495
496 CALL timeset(routinen, handle)
497 NULLIFY (scon_oba, scon_obb, scon_fba, scon_fbb, scona_mix, sconb_mix, &
498 oba_index, obb_index, fba_index, fbb_index)
499 dab = sqrt(sum(rab*rab))
500 nba = obasa%nsgf
501 nbb = obasb%nsgf
502 nfa = fbasa%nsgf
503 nfb = fbasb%nsgf
504
505 !*** calculate overlap integrals; for iatom=jatom this is the self-overlap
506 IF (iatom == jatom .AND. dab < lri_env%delta) THEN
507 !*** (fa,fa)
508 lrii%sab(1:nfa, 1:nfa) = lri_env%bas_prop(ikind)%ri_ovlp(1:nfa, 1:nfa)
509 !*** (a,a)
510 lrii%soo(1:nba, 1:nba) = lri_env%bas_prop(ikind)%orb_ovlp(1:nba, 1:nba)
511 !*** (a,a,fa)
512 lrii%abaint(1:nba, 1:nba, 1:nfa) = lri_env%bas_prop(ikind)%ovlp3
513 IF (calculate_forces) THEN
514 lrii%dsab = 0._dp
515 lrii%dsoo = 0._dp
516 lrii%dabdaint = 0.0_dp
517 lrii%dabbint = 0.0_dp
518 END IF
519 ELSE
520 IF (calculate_forces) THEN
521 do_force = lrii%calc_force_pair
522 ELSE
523 do_force = .false.
524 END IF
525 scon_oba => lri_env%bas_prop(ikind)%scon_orb
526 scon_obb => lri_env%bas_prop(jkind)%scon_orb
527 scon_fba => lri_env%bas_prop(ikind)%scon_ri
528 scon_fbb => lri_env%bas_prop(jkind)%scon_ri
529 scona_mix => lri_env%bas_prop(ikind)%scon_mix
530 sconb_mix => lri_env%bas_prop(jkind)%scon_mix
531 oba_index => lri_env%bas_prop(ikind)%orb_index
532 fba_index => lri_env%bas_prop(ikind)%ri_index
533 obb_index => lri_env%bas_prop(jkind)%orb_index
534 fbb_index => lri_env%bas_prop(jkind)%ri_index
535 CALL lri_precalc_angular_shg_part(obasa, obasb, fbasa, fbasb, rab, waux_mat, dwaux_mat, &
536 do_force)
537 !*** (fa,fb)
538 CALL int_overlap_ab_shg_low(lrii%sab, lrii%dsab, rab, fbasa, fbasb, scon_fba, scon_fbb, &
539 waux_mat, dwaux_mat, .true., do_force, contraction_high=.false.)
540 !*** (a,b)
541 CALL int_overlap_ab_shg_low(lrii%soo, lrii%dsoo, rab, obasa, obasb, scon_oba, scon_obb, &
542 waux_mat, dwaux_mat, .true., do_force, contraction_high=.true.)
543 !*** (a,b,fa)
544 CALL int_overlap_aba_shg_low(lrii%abaint, lrii%dabdaint, rab, obasa, obasb, fbasa, &
545 scon_obb, scona_mix, oba_index, fba_index, &
546 lri_env%cg_shg%cg_coeff, lri_env%cg_shg%cg_none0_list, &
547 lri_env%cg_shg%ncg_none0, &
548 waux_mat, dwaux_mat, .true., do_force)
549 !*** (a,b,fb)
550 IF (ikind == jkind) THEN
551 CALL get_abb_same_kind(lrii%abbint, lrii%dabbint, lrii%abaint, lrii%dabdaint, &
552 rab, obasa, fbasa, .true., do_force)
553 ELSE
554 CALL int_overlap_abb_shg_low(lrii%abbint, lrii%dabbint, rab, obasa, obasb, fbasb, &
555 scon_oba, sconb_mix, obb_index, fbb_index, &
556 lri_env%cg_shg%cg_coeff, lri_env%cg_shg%cg_none0_list, &
557 lri_env%cg_shg%ncg_none0, &
558 waux_mat, dwaux_mat, .true., do_force)
559 END IF
560
561 DEALLOCATE (waux_mat, dwaux_mat)
562 END IF
563
564 CALL timestop(handle)
565
566 END SUBROUTINE lri_int_shg
567
568! **************************************************************************************************
569!> \brief ...
570!> \param lriint ...
571!> \param lridint ...
572!> \param nba ...
573!> \param nbb ...
574!> \param nfa ...
575!> \param nfb ...
576!> \param skip_sab ...
577!> \param skip_soo ...
578!> \param skip_aba ...
579!> \param skip_abb ...
580!> \param skip_dsab ...
581!> \param skip_dsoo ...
582!> \param skip_daba ...
583!> \param skip_dabb ...
584! **************************************************************************************************
585 SUBROUTINE allocate_int_type(lriint, lridint, nba, nbb, nfa, nfb, &
586 skip_sab, skip_soo, skip_aba, skip_abb, &
587 skip_dsab, skip_dsoo, skip_daba, skip_dabb)
588
589 TYPE(int_type), INTENT(INOUT), OPTIONAL :: lriint
590 TYPE(dint_type), INTENT(INOUT), OPTIONAL :: lridint
591 INTEGER, INTENT(IN) :: nba, nbb, nfa, nfb
592 LOGICAL, INTENT(IN), OPTIONAL :: skip_sab, skip_soo, skip_aba, skip_abb, &
593 skip_dsab, skip_dsoo, skip_daba, &
594 skip_dabb
595
596 LOGICAL :: do_aba, do_abb, do_daba, do_dabb, &
597 do_dsab, do_dsoo, do_sab, do_soo
598
599 IF (PRESENT(lriint)) THEN
600 do_sab = .true.
601 IF (PRESENT(skip_sab)) do_sab = .NOT. skip_sab
602 do_soo = .true.
603 IF (PRESENT(skip_soo)) do_soo = .NOT. skip_soo
604 do_aba = .true.
605 IF (PRESENT(skip_aba)) do_aba = .NOT. skip_aba
606 do_abb = .true.
607 IF (PRESENT(skip_abb)) do_abb = .NOT. skip_abb
608 !
609 IF (do_sab) THEN
610 IF (ALLOCATED(lriint%sabint)) DEALLOCATE (lriint%sabint)
611 ALLOCATE (lriint%sabint(nfa, nfb))
612 lriint%sabint = 0.0_dp
613 END IF
614 IF (do_soo) THEN
615 IF (ALLOCATED(lriint%sooint)) DEALLOCATE (lriint%sooint)
616 ALLOCATE (lriint%sooint(nba, nbb))
617 lriint%sooint = 0.0_dp
618 END IF
619 IF (do_aba) THEN
620 IF (ALLOCATED(lriint%abaint)) DEALLOCATE (lriint%abaint)
621 ALLOCATE (lriint%abaint(nba, nbb, nfa))
622 lriint%abaint = 0.0_dp
623 END IF
624 IF (do_abb) THEN
625 IF (ALLOCATED(lriint%abbint)) DEALLOCATE (lriint%abbint)
626 ALLOCATE (lriint%abbint(nba, nbb, nfb))
627 lriint%abbint = 0.0_dp
628 END IF
629 END IF
630 !
631 IF (PRESENT(lridint)) THEN
632 do_dsab = .true.
633 IF (PRESENT(skip_dsab)) do_dsab = .NOT. skip_dsab
634 do_dsoo = .true.
635 IF (PRESENT(skip_dsoo)) do_dsoo = .NOT. skip_dsoo
636 do_daba = .true.
637 IF (PRESENT(skip_daba)) do_daba = .NOT. skip_daba
638 do_dabb = .true.
639 IF (PRESENT(skip_dabb)) do_dabb = .NOT. skip_dabb
640 !
641 IF (do_dsab) THEN
642 IF (ALLOCATED(lridint%dsabint)) DEALLOCATE (lridint%dsabint)
643 ALLOCATE (lridint%dsabint(nfa, nfb, 3))
644 lridint%dsabint = 0.0_dp
645 END IF
646 IF (do_dsoo) THEN
647 IF (ALLOCATED(lridint%dsooint)) DEALLOCATE (lridint%dsooint)
648 ALLOCATE (lridint%dsooint(nba, nbb, 3))
649 lridint%dsooint = 0.0_dp
650 END IF
651 IF (do_daba) THEN
652 IF (ALLOCATED(lridint%dabaint)) DEALLOCATE (lridint%dabaint)
653 ALLOCATE (lridint%dabaint(nba, nbb, nfa, 3))
654 lridint%dabaint = 0.0_dp
655 END IF
656 IF (do_dabb) THEN
657 IF (ALLOCATED(lridint%dabbint)) DEALLOCATE (lridint%dabbint)
658 ALLOCATE (lridint%dabbint(nba, nbb, nfb, 3))
659 lridint%dabbint = 0.0_dp
660 END IF
661 END IF
662
663 END SUBROUTINE allocate_int_type
664
665! **************************************************************************************************
666!> \brief ...
667!> \param lriint ...
668!> \param lridint ...
669! **************************************************************************************************
670 SUBROUTINE deallocate_int_type(lriint, lridint)
671
672 TYPE(int_type), INTENT(INOUT), OPTIONAL :: lriint
673 TYPE(dint_type), INTENT(INOUT), OPTIONAL :: lridint
674
675 IF (PRESENT(lriint)) THEN
676 IF (ALLOCATED(lriint%sabint)) DEALLOCATE (lriint%sabint)
677 IF (ALLOCATED(lriint%sooint)) DEALLOCATE (lriint%sooint)
678 IF (ALLOCATED(lriint%abaint)) DEALLOCATE (lriint%abaint)
679 IF (ALLOCATED(lriint%abbint)) DEALLOCATE (lriint%abbint)
680 END IF
681 !
682 IF (PRESENT(lridint)) THEN
683 IF (ALLOCATED(lridint%dsabint)) DEALLOCATE (lridint%dsabint)
684 IF (ALLOCATED(lridint%dsooint)) DEALLOCATE (lridint%dsooint)
685 IF (ALLOCATED(lridint%dabaint)) DEALLOCATE (lridint%dabaint)
686 IF (ALLOCATED(lridint%dabbint)) DEALLOCATE (lridint%dabbint)
687 END IF
688
689 END SUBROUTINE deallocate_int_type
690
691END MODULE lri_integrals
basis_set_types
Definition basis_set_types.F:15
generic_os_integrals
Calculation of contracted, spherical Gaussian integrals using the (OS) integral scheme....
Definition generic_os_integrals.F:18
generic_os_integrals::int_overlap_aba_os
subroutine, public int_overlap_aba_os(abaint, dabdaint, rab, oba, obb, fba, calculate_forces, debug, dmax)
calculate integrals (a,b,fa)
Definition generic_os_integrals.F:506
generic_os_integrals::int_overlap_ab_os
subroutine, public int_overlap_ab_os(sab, dsab, rab, fba, fbb, calculate_forces, debug, dmax)
calculate overlap integrals (a,b)
Definition generic_os_integrals.F:270
generic_os_integrals::int_overlap_abb_os
subroutine, public int_overlap_abb_os(abbint, dabbint, rab, oba, obb, fbb, calculate_forces, debug, dmax)
calculate integrals (a,b,fb)
Definition generic_os_integrals.F:667
generic_shg_integrals
Calculation of contracted, spherical Gaussian integrals using the solid harmonic Gaussian (SHG) integ...
Definition generic_shg_integrals.F:21
generic_shg_integrals::lri_precalc_angular_shg_part
subroutine, public lri_precalc_angular_shg_part(oba, obb, fba, fbb, rab, waux_mat, dwaux_mat, calculate_forces)
precalculates the angular part of the SHG integrals for the matrices (fa,fb), (a,b),...
Definition generic_shg_integrals.F:1242
generic_shg_integrals::int_overlap_aba_shg_low
subroutine, public int_overlap_aba_shg_low(abaint, dabdaint, rab, oba, obb, fba, scon_obb, scona_mix, oba_index, fba_index, cg_coeff, cg_none0_list, ncg_none0, waux_mat, dwaux_mat, calculate_ints, calculate_forces)
calculate integrals (a,b,fa); requires angular-dependent part as input
Definition generic_shg_integrals.F:1078
generic_shg_integrals::int_overlap_abb_shg_low
subroutine, public int_overlap_abb_shg_low(abbint, dabbint, rab, oba, obb, fbb, scon_oba, sconb_mix, obb_index, fbb_index, cg_coeff, cg_none0_list, ncg_none0, waux_mat, dwaux_mat, calculate_ints, calculate_forces)
calculate integrals (a,b,fb); requires angular-dependent part as input
Definition generic_shg_integrals.F:772
generic_shg_integrals::get_abb_same_kind
subroutine, public get_abb_same_kind(abbint, dabbint, abaint, dabdaint, rab, oba, fba, calculate_ints, calculate_forces)
obtain integrals (a,b,fb) by symmetry relations from (a,b,fa) if basis sets at a and b are of the sam...
Definition generic_shg_integrals.F:935
generic_shg_integrals::int_overlap_ab_shg_low
subroutine, public int_overlap_ab_shg_low(sab, dsab, rab, fba, fbb, scona_shg, sconb_shg, waux_mat, dwaux_mat, calculate_ints, calculate_forces, contraction_high)
calculate overlap integrals (a,b); requires angular-dependent part as input
Definition generic_shg_integrals.F:526
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
lri_environment_types
contains the types and subroutines for dealing with the lri_env lri : local resolution of the identit...
Definition lri_environment_types.F:18
lri_integrals
Calculates integrals for LRIGPW method lri : local resolution of the identity.
Definition lri_integrals.F:19
lri_integrals::deallocate_int_type
subroutine, public deallocate_int_type(lriint, lridint)
...
Definition lri_integrals.F:671
lri_integrals::allocate_int_type
subroutine, public allocate_int_type(lriint, lridint, nba, nbb, nfa, nfb, skip_sab, skip_soo, skip_aba, skip_abb, skip_dsab, skip_dsoo, skip_daba, skip_dabb)
...
Definition lri_integrals.F:588
lri_integrals::lri_dint
subroutine, public lri_dint(lri_env, lrii, rab, obasa, obasb, fbasa, fbasb, iatom, jatom, ikind, jkind, calculate_forces)
...
Definition lri_integrals.F:233
lri_integrals::lri_int2
subroutine, public lri_int2(lri_env, lrii, lriint, rab, obasa, obasb, fbasa, fbasb, iatom, jatom, ikind, jkind)
calcuates the lri integrals using solid harmonic Gaussians
Definition lri_integrals.F:113
lri_integrals::lri_int
subroutine, public lri_int(lri_env, lrii, rab, obasa, obasb, fbasa, fbasb, iatom, jatom, ikind, jkind, calculate_forces)
calcuates the lri integrals using solid harmonic Gaussians
Definition lri_integrals.F:79
lri_integrals::lri_dint2
subroutine, public lri_dint2(lri_env, lrii, lridint, rab, obasa, obasb, fbasa, fbasb, iatom, jatom, ikind, jkind)
...
Definition lri_integrals.F:268
basis_set_types::gto_basis_set_type
Definition basis_set_types.F:71
lri_environment_types::lri_environment_type
Definition lri_environment_types.F:269
lri_environment_types::lri_int_type
Definition lri_environment_types.F:91
lri_integrals::dint_type
Definition lri_integrals.F:46
lri_integrals::int_type
Definition lri_integrals.F:39