Line data Source code
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 : !> \brief 3-center overlap type integrals containers
9 : !> \par History
10 : !> - Added options to only keep (abc) triplet if b and c share the same center (2019 A.Bussy)
11 : ! **************************************************************************************************
12 : MODULE qs_o3c_types
13 :
14 : USE basis_set_types, ONLY: gto_basis_set_p_type
15 : USE kinds, ONLY: dp
16 : USE qs_neighbor_list_types, ONLY: &
17 : get_iterator_info, get_neighbor_list_set_p, neighbor_list_iterate, &
18 : neighbor_list_iterator_create, neighbor_list_iterator_p_type, &
19 : neighbor_list_iterator_release, neighbor_list_set_p_type, nl_set_sub_iterator, &
20 : nl_sub_iterate
21 : #include "./base/base_uses.f90"
22 :
23 : IMPLICIT NONE
24 :
25 : PRIVATE
26 :
27 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_o3c_types'
28 :
29 : ! **************************************************************************************************
30 : ! O3C Integrals
31 : ! **************************************************************************************************
32 :
33 : TYPE o3c_int_type
34 : PRIVATE
35 : INTEGER :: katom = -1, kkind = -1
36 : INTEGER :: ni = -1, nj = -1, nk = -1
37 : REAL(KIND=dp), DIMENSION(3) :: rik = -1.0_dp
38 : INTEGER, DIMENSION(3) :: cellk = -1
39 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: integral => NULL()
40 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: tvec => NULL()
41 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: force_i => NULL()
42 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: force_j => NULL()
43 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: force_k => NULL()
44 : END TYPE o3c_int_type
45 :
46 : TYPE o3c_pair_type
47 : PRIVATE
48 : INTEGER :: iatom = -1, ikind = -1
49 : INTEGER :: jatom = -1, jkind = -1
50 : REAL(KIND=dp), DIMENSION(3) :: rij = -1.0_dp
51 : INTEGER, DIMENSION(3) :: cellj = -1
52 : INTEGER :: nklist = -1
53 : TYPE(o3c_int_type), DIMENSION(:), POINTER :: ijk => NULL()
54 : END TYPE o3c_pair_type
55 :
56 : TYPE o3c_container_type
57 : PRIVATE
58 : LOGICAL :: ijsymmetric = .FALSE.
59 : INTEGER :: nijpairs = -1
60 : INTEGER :: nspin = -1
61 : TYPE(o3c_pair_type), DIMENSION(:), POINTER :: ijpair => NULL()
62 : ! basis sets and neighbor lists are pointing to other resources
63 : ! we don't keep track if the data is available and correct
64 : TYPE(gto_basis_set_p_type), DIMENSION(:), &
65 : POINTER :: basis_set_list_a => NULL(), basis_set_list_b => NULL(), &
66 : basis_set_list_c => NULL()
67 : TYPE(neighbor_list_set_p_type), &
68 : DIMENSION(:), POINTER :: sab_nl => NULL(), sac_nl => NULL()
69 : END TYPE o3c_container_type
70 :
71 : ! **************************************************************************************************
72 : ! O3C Iterator
73 : ! **************************************************************************************************
74 :
75 : TYPE o3c_iterator_type
76 : PRIVATE
77 : TYPE(o3c_container_type), POINTER :: o3c => NULL()
78 : INTEGER :: ijp_last = -1, k_last = -1
79 : INTEGER, DIMENSION(:), POINTER :: ijp_thread => NULL(), k_thread => NULL()
80 : END TYPE o3c_iterator_type
81 :
82 : ! **************************************************************************************************
83 : ! O3C vector
84 : ! **************************************************************************************************
85 :
86 : TYPE o3c_vec_type
87 : PRIVATE
88 : INTEGER :: n = -1
89 : REAL(KIND=dp), DIMENSION(:), POINTER :: v => NULL()
90 : END TYPE o3c_vec_type
91 :
92 : ! **************************************************************************************************
93 :
94 : PUBLIC :: o3c_container_type
95 : PUBLIC :: release_o3c_container, init_o3c_container, get_o3c_container, set_o3c_container
96 : PUBLIC :: o3c_iterator_type
97 : PUBLIC :: o3c_iterator_create, o3c_iterator_release, get_o3c_iterator_info, o3c_iterate
98 : PUBLIC :: o3c_vec_type
99 : PUBLIC :: o3c_vec_create, o3c_vec_release, get_o3c_vec
100 :
101 : CONTAINS
102 :
103 : ! **************************************************************************************************
104 : !> \brief ...
105 : !> \param o3c ...
106 : !> \param nspin ...
107 : !> \param basis_set_list_a ...
108 : !> \param basis_set_list_b ...
109 : !> \param basis_set_list_c ...
110 : !> \param sab_nl ...
111 : !> \param sac_nl ...
112 : !> \param only_bc_same_center only consider a,b,c atoms if b and c share the same center
113 : !> \par History: only_bc_same_cetner added by A.Bussy for XAS_TDP (04.2019)
114 : ! **************************************************************************************************
115 132 : SUBROUTINE init_o3c_container(o3c, nspin, basis_set_list_a, basis_set_list_b, basis_set_list_c, &
116 : sab_nl, sac_nl, only_bc_same_center)
117 : TYPE(o3c_container_type) :: o3c
118 : INTEGER, INTENT(IN) :: nspin
119 : TYPE(gto_basis_set_p_type), DIMENSION(:), POINTER :: basis_set_list_a, basis_set_list_b, &
120 : basis_set_list_c
121 : TYPE(neighbor_list_set_p_type), DIMENSION(:), &
122 : POINTER :: sab_nl, sac_nl
123 : LOGICAL, INTENT(IN), OPTIONAL :: only_bc_same_center
124 :
125 : INTEGER :: kkind, nij, nk, nkind
126 : LOGICAL :: my_sort_bc, symmetric
127 : REAL(dp) :: rik(3), rjk(3)
128 : TYPE(neighbor_list_iterator_p_type), &
129 132 : DIMENSION(:), POINTER :: ac_iterator, nl_iterator
130 : TYPE(o3c_int_type), POINTER :: ijk
131 : TYPE(o3c_pair_type), POINTER :: ijpair
132 :
133 132 : CALL get_neighbor_list_set_p(sab_nl, symmetric=symmetric)
134 132 : o3c%ijsymmetric = symmetric
135 132 : CPASSERT(symmetric)
136 :
137 132 : o3c%nspin = nspin
138 :
139 132 : o3c%basis_set_list_a => basis_set_list_a
140 132 : o3c%basis_set_list_b => basis_set_list_b
141 132 : o3c%basis_set_list_c => basis_set_list_c
142 :
143 132 : o3c%sab_nl => sab_nl
144 132 : o3c%sac_nl => sac_nl
145 :
146 132 : nkind = SIZE(basis_set_list_a)
147 :
148 132 : my_sort_bc = .FALSE.
149 132 : IF (PRESENT(only_bc_same_center)) my_sort_bc = only_bc_same_center
150 :
151 : ! determine the number of ij pairs
152 132 : nij = 0
153 132 : CALL neighbor_list_iterator_create(nl_iterator, sab_nl)
154 17360 : DO WHILE (neighbor_list_iterate(nl_iterator) == 0)
155 17228 : nij = nij + 1
156 : END DO
157 132 : CALL neighbor_list_iterator_release(nl_iterator)
158 132 : o3c%nijpairs = nij
159 132 : NULLIFY (o3c%ijpair)
160 18528 : ALLOCATE (o3c%ijpair(nij))
161 :
162 : ! for each pair set up the ijk lists
163 132 : nij = 0
164 132 : CALL neighbor_list_iterator_create(nl_iterator, sab_nl)
165 132 : CALL neighbor_list_iterator_create(ac_iterator, sac_nl, search=.TRUE.)
166 17360 : DO WHILE (neighbor_list_iterate(nl_iterator) == 0)
167 17228 : nij = nij + 1
168 17228 : ijpair => o3c%ijpair(nij)
169 : CALL get_iterator_info(nl_iterator, ikind=ijpair%ikind, jkind=ijpair%jkind, &
170 : iatom=ijpair%iatom, jatom=ijpair%jatom, &
171 17228 : r=ijpair%rij, cell=ijpair%cellj)
172 17228 : NULLIFY (ijpair%ijk)
173 17228 : nk = 0
174 50446 : DO kkind = 1, nkind
175 33218 : CALL nl_set_sub_iterator(ac_iterator, ijpair%ikind, kkind, ijpair%iatom)
176 71754 : DO WHILE (nl_sub_iterate(ac_iterator) == 0)
177 21308 : IF (my_sort_bc) THEN
178 : !we only take ijk if rjk = 0 OR rik = 0 (because of symmetry)
179 690 : CALL get_iterator_info(ac_iterator, r=rik)
180 2760 : rjk(:) = rik(:) - ijpair%rij(:)
181 2979 : IF (.NOT. (ALL(ABS(rjk) .LE. 1.0E-4_dp) .OR. ALL(ABS(rik) .LE. 1.0E-4_dp))) CYCLE
182 : END IF
183 21308 : nk = nk + 1
184 : END DO
185 : END DO
186 : ! ijk lists
187 17228 : ijpair%nklist = nk
188 190249 : ALLOCATE (ijpair%ijk(nk))
189 : ! fill the ijk lists
190 17228 : nk = 0
191 50578 : DO kkind = 1, nkind
192 33218 : CALL nl_set_sub_iterator(ac_iterator, ijpair%ikind, kkind, ijpair%iatom)
193 71754 : DO WHILE (nl_sub_iterate(ac_iterator) == 0)
194 21308 : IF (my_sort_bc) THEN
195 : !we only take ijk if rjk = 0 OR rik = 0 (because of symmetry)
196 690 : CALL get_iterator_info(ac_iterator, r=rik)
197 2760 : rjk(:) = rik(:) - ijpair%rij(:)
198 2979 : IF (.NOT. (ALL(ABS(rjk) .LE. 1.0E-4_dp) .OR. ALL(ABS(rik) .LE. 1.0E-4_dp))) CYCLE
199 : END IF
200 :
201 21018 : nk = nk + 1
202 21018 : ijk => ijpair%ijk(nk)
203 21018 : CALL get_iterator_info(ac_iterator, jatom=ijk%katom, r=ijk%rik, cell=ijk%cellk)
204 21018 : ijk%kkind = kkind
205 21018 : ijk%ni = 0
206 21018 : ijk%nj = 0
207 21018 : ijk%nk = 0
208 21018 : NULLIFY (ijk%integral)
209 21018 : NULLIFY (ijk%tvec)
210 21018 : NULLIFY (ijk%force_i)
211 21018 : NULLIFY (ijk%force_j)
212 21308 : NULLIFY (ijk%force_k)
213 : END DO
214 : END DO
215 : END DO
216 132 : CALL neighbor_list_iterator_release(ac_iterator)
217 132 : CALL neighbor_list_iterator_release(nl_iterator)
218 :
219 132 : END SUBROUTINE init_o3c_container
220 : ! **************************************************************************************************
221 : !> \brief ...
222 : !> \param o3c_container ...
223 : ! **************************************************************************************************
224 132 : SUBROUTINE release_o3c_container(o3c_container)
225 :
226 : TYPE(o3c_container_type) :: o3c_container
227 :
228 132 : o3c_container%ijsymmetric = .FALSE.
229 132 : o3c_container%nijpairs = 0
230 :
231 132 : NULLIFY (o3c_container%basis_set_list_a)
232 132 : NULLIFY (o3c_container%basis_set_list_b)
233 132 : NULLIFY (o3c_container%basis_set_list_c)
234 :
235 132 : NULLIFY (o3c_container%sab_nl)
236 132 : NULLIFY (o3c_container%sac_nl)
237 :
238 132 : IF (ASSOCIATED(o3c_container%ijpair)) THEN
239 132 : CALL release_ijpair(o3c_container%ijpair)
240 132 : DEALLOCATE (o3c_container%ijpair)
241 : END IF
242 :
243 132 : END SUBROUTINE release_o3c_container
244 :
245 : ! **************************************************************************************************
246 : !> \brief ...
247 : !> \param ijpair ...
248 : ! **************************************************************************************************
249 132 : SUBROUTINE release_ijpair(ijpair)
250 :
251 : TYPE(o3c_pair_type), DIMENSION(:) :: ijpair
252 :
253 : INTEGER :: i
254 :
255 17360 : DO i = 1, SIZE(ijpair)
256 17228 : ijpair(i)%iatom = 0
257 17228 : ijpair(i)%ikind = 0
258 17228 : ijpair(i)%jatom = 0
259 17228 : ijpair(i)%jkind = 0
260 17228 : ijpair(i)%nklist = 0
261 68912 : ijpair(i)%rij = 0.0_dp
262 68912 : ijpair(i)%cellj = 0
263 17360 : IF (ASSOCIATED(ijpair(i)%ijk)) THEN
264 17228 : CALL release_ijk(ijpair(i)%ijk)
265 17228 : DEALLOCATE (ijpair(i)%ijk)
266 : END IF
267 : END DO
268 :
269 132 : END SUBROUTINE release_ijpair
270 :
271 : ! **************************************************************************************************
272 : !> \brief ...
273 : !> \param ijk ...
274 : ! **************************************************************************************************
275 17228 : SUBROUTINE release_ijk(ijk)
276 :
277 : TYPE(o3c_int_type), DIMENSION(:) :: ijk
278 :
279 : INTEGER :: i
280 :
281 38246 : DO i = 1, SIZE(ijk)
282 21018 : ijk(i)%katom = 0
283 21018 : ijk(i)%kkind = 0
284 21018 : ijk(i)%ni = 0
285 21018 : ijk(i)%nj = 0
286 21018 : ijk(i)%nk = 0
287 84072 : ijk(i)%rik = 0.0_dp
288 84072 : ijk(i)%cellk = 0
289 21018 : IF (ASSOCIATED(ijk(i)%integral)) THEN
290 0 : DEALLOCATE (ijk(i)%integral)
291 : END IF
292 21018 : IF (ASSOCIATED(ijk(i)%tvec)) THEN
293 0 : DEALLOCATE (ijk(i)%tvec)
294 : END IF
295 21018 : IF (ASSOCIATED(ijk(i)%force_i)) THEN
296 0 : DEALLOCATE (ijk(i)%force_i)
297 : END IF
298 21018 : IF (ASSOCIATED(ijk(i)%force_j)) THEN
299 0 : DEALLOCATE (ijk(i)%force_j)
300 : END IF
301 38246 : IF (ASSOCIATED(ijk(i)%force_k)) THEN
302 0 : DEALLOCATE (ijk(i)%force_k)
303 : END IF
304 : END DO
305 :
306 17228 : END SUBROUTINE release_ijk
307 :
308 : ! **************************************************************************************************
309 : !> \brief ...
310 : !> \param o3c ...
311 : !> \param ijsymmetric ...
312 : !> \param nspin ...
313 : !> \param nijpairs ...
314 : !> \param ijpair ...
315 : !> \param basis_set_list_a ...
316 : !> \param basis_set_list_b ...
317 : !> \param basis_set_list_c ...
318 : !> \param sab_nl ...
319 : !> \param sac_nl ...
320 : ! **************************************************************************************************
321 :
322 0 : SUBROUTINE get_o3c_container(o3c, ijsymmetric, nspin, nijpairs, ijpair, &
323 : basis_set_list_a, basis_set_list_b, basis_set_list_c, &
324 : sab_nl, sac_nl)
325 : TYPE(o3c_container_type) :: o3c
326 : LOGICAL, OPTIONAL :: ijsymmetric
327 : INTEGER, OPTIONAL :: nspin, nijpairs
328 : TYPE(o3c_pair_type), DIMENSION(:), OPTIONAL, &
329 : POINTER :: ijpair
330 : TYPE(gto_basis_set_p_type), DIMENSION(:), &
331 : OPTIONAL, POINTER :: basis_set_list_a, basis_set_list_b, &
332 : basis_set_list_c
333 : TYPE(neighbor_list_set_p_type), DIMENSION(:), &
334 : OPTIONAL, POINTER :: sab_nl, sac_nl
335 :
336 0 : IF (PRESENT(ijsymmetric)) ijsymmetric = o3c%ijsymmetric
337 0 : IF (PRESENT(nspin)) nspin = o3c%nspin
338 0 : IF (PRESENT(nijpairs)) nijpairs = o3c%nijpairs
339 0 : IF (PRESENT(ijpair)) ijpair => o3c%ijpair
340 0 : IF (PRESENT(basis_set_list_a)) basis_set_list_a => o3c%basis_set_list_a
341 0 : IF (PRESENT(basis_set_list_b)) basis_set_list_b => o3c%basis_set_list_b
342 0 : IF (PRESENT(basis_set_list_c)) basis_set_list_c => o3c%basis_set_list_c
343 0 : IF (PRESENT(sab_nl)) sab_nl => o3c%sab_nl
344 0 : IF (PRESENT(sac_nl)) sac_nl => o3c%sac_nl
345 :
346 0 : END SUBROUTINE get_o3c_container
347 :
348 : ! **************************************************************************************************
349 : ! O3C Iterator
350 : ! **************************************************************************************************
351 : !> \brief ...
352 : !> \param o3c ...
353 : !> \param o3c_iterator ...
354 : !> \param nthread ...
355 : ! **************************************************************************************************
356 132 : SUBROUTINE o3c_iterator_create(o3c, o3c_iterator, nthread)
357 : TYPE(o3c_container_type), POINTER :: o3c
358 : TYPE(o3c_iterator_type) :: o3c_iterator
359 : INTEGER, OPTIONAL :: nthread
360 :
361 : INTEGER :: n
362 :
363 132 : IF (PRESENT(nthread)) THEN
364 132 : n = nthread
365 : ELSE
366 : n = 1
367 : END IF
368 :
369 132 : o3c_iterator%o3c => o3c
370 132 : o3c_iterator%ijp_last = 0
371 132 : o3c_iterator%k_last = 0
372 396 : ALLOCATE (o3c_iterator%ijp_thread(0:n - 1))
373 396 : ALLOCATE (o3c_iterator%k_thread(0:n - 1))
374 264 : o3c_iterator%ijp_thread = 0
375 264 : o3c_iterator%k_thread = 0
376 :
377 132 : END SUBROUTINE o3c_iterator_create
378 :
379 : ! **************************************************************************************************
380 : !> \brief ...
381 : !> \param o3c_iterator ...
382 : ! **************************************************************************************************
383 132 : SUBROUTINE o3c_iterator_release(o3c_iterator)
384 : TYPE(o3c_iterator_type) :: o3c_iterator
385 :
386 132 : NULLIFY (o3c_iterator%o3c)
387 132 : o3c_iterator%ijp_last = 0
388 132 : o3c_iterator%k_last = 0
389 132 : DEALLOCATE (o3c_iterator%ijp_thread)
390 132 : DEALLOCATE (o3c_iterator%k_thread)
391 :
392 132 : END SUBROUTINE o3c_iterator_release
393 :
394 : ! **************************************************************************************************
395 : !> \brief ...
396 : !> \param o3c_iterator ...
397 : !> \param mepos ...
398 : !> \param iatom ...
399 : !> \param jatom ...
400 : !> \param katom ...
401 : !> \param ikind ...
402 : !> \param jkind ...
403 : !> \param kkind ...
404 : !> \param rij ...
405 : !> \param rik ...
406 : !> \param cellj ...
407 : !> \param cellk ...
408 : !> \param integral ...
409 : !> \param tvec ...
410 : !> \param force_i ...
411 : !> \param force_j ...
412 : !> \param force_k ...
413 : ! **************************************************************************************************
414 21018 : SUBROUTINE get_o3c_iterator_info(o3c_iterator, mepos, &
415 : iatom, jatom, katom, ikind, jkind, kkind, &
416 : rij, rik, cellj, cellk, &
417 : integral, tvec, force_i, force_j, force_k)
418 : TYPE(o3c_iterator_type) :: o3c_iterator
419 : INTEGER, OPTIONAL :: mepos, iatom, jatom, katom, ikind, &
420 : jkind, kkind
421 : REAL(KIND=dp), DIMENSION(3), OPTIONAL :: rij, rik
422 : INTEGER, DIMENSION(3), OPTIONAL :: cellj, cellk
423 : REAL(KIND=dp), DIMENSION(:, :, :), OPTIONAL, &
424 : POINTER :: integral
425 : REAL(KIND=dp), DIMENSION(:, :), OPTIONAL, POINTER :: tvec, force_i, force_j, force_k
426 :
427 : INTEGER :: ij, k, me
428 : TYPE(o3c_container_type), POINTER :: o3c
429 : TYPE(o3c_int_type), POINTER :: ijk
430 : TYPE(o3c_pair_type), POINTER :: ijp
431 :
432 21018 : IF (PRESENT(mepos)) THEN
433 21018 : me = mepos
434 : ELSE
435 : me = 0
436 : END IF
437 :
438 21018 : ij = o3c_iterator%ijp_thread(me)
439 21018 : k = o3c_iterator%k_thread(me)
440 :
441 21018 : o3c => o3c_iterator%o3c
442 21018 : ijp => o3c%ijpair(ij)
443 21018 : ijk => ijp%ijk(k)
444 :
445 21018 : IF (PRESENT(iatom)) iatom = ijp%iatom
446 21018 : IF (PRESENT(jatom)) jatom = ijp%jatom
447 21018 : IF (PRESENT(ikind)) ikind = ijp%ikind
448 21018 : IF (PRESENT(jkind)) jkind = ijp%jkind
449 21018 : IF (PRESENT(katom)) katom = ijk%katom
450 21018 : IF (PRESENT(kkind)) kkind = ijk%kkind
451 :
452 105090 : IF (PRESENT(rij)) rij(1:3) = ijp%rij(1:3)
453 105090 : IF (PRESENT(rik)) rik(1:3) = ijk%rik(1:3)
454 :
455 21018 : IF (PRESENT(cellj)) cellj(1:3) = ijp%cellj(1:3)
456 21018 : IF (PRESENT(cellk)) cellk(1:3) = ijk%cellk(1:3)
457 :
458 21018 : IF (PRESENT(integral)) integral => ijk%integral
459 21018 : IF (PRESENT(tvec)) tvec => ijk%tvec
460 21018 : IF (PRESENT(force_i)) force_i => ijk%force_i
461 21018 : IF (PRESENT(force_j)) force_j => ijk%force_j
462 21018 : IF (PRESENT(force_k)) force_k => ijk%force_k
463 :
464 21018 : END SUBROUTINE get_o3c_iterator_info
465 :
466 : ! **************************************************************************************************
467 : !> \brief ...
468 : !> \param o3c_iterator ...
469 : !> \param mepos ...
470 : !> \return ...
471 : ! **************************************************************************************************
472 21150 : FUNCTION o3c_iterate(o3c_iterator, mepos) RESULT(istat)
473 : TYPE(o3c_iterator_type) :: o3c_iterator
474 : INTEGER, OPTIONAL :: mepos
475 : INTEGER :: istat
476 :
477 : INTEGER :: ij, ijpair, klist, me
478 : TYPE(o3c_container_type), POINTER :: o3c
479 :
480 21150 : IF (PRESENT(mepos)) THEN
481 21150 : me = mepos
482 : ELSE
483 : me = 0
484 : END IF
485 :
486 : !If the neighbors lists are restricted (XAS_TDP), might have nijpairs = 0 on some procs
487 21150 : IF (o3c_iterator%o3c%nijpairs == 0) THEN
488 21150 : istat = 1
489 : RETURN
490 : END IF
491 :
492 42260 : !$OMP CRITICAL(o3c_iterate_critical)
493 21130 : o3c => o3c_iterator%o3c
494 : ! we iterate from the last position
495 21130 : ijpair = o3c_iterator%ijp_last
496 21130 : klist = o3c_iterator%k_last
497 :
498 21130 : IF (ijpair == 0 .AND. klist == 0) THEN
499 : ! first step
500 112 : istat = 1
501 131 : DO ij = 1, o3c%nijpairs
502 131 : IF (o3c%ijpair(ij)%nklist > 0) THEN
503 112 : o3c_iterator%ijp_thread(me) = ij
504 112 : o3c_iterator%k_thread(me) = 1
505 : istat = 0
506 : EXIT
507 : END IF
508 : END DO
509 21018 : ELSE IF (ijpair == o3c%nijpairs .AND. klist == o3c%ijpair(ijpair)%nklist) THEN
510 : ! last step reached
511 : istat = 1
512 20915 : ELSE IF (klist == o3c%ijpair(ijpair)%nklist) THEN
513 : ! last step in this ij list
514 14076 : istat = 1
515 17106 : DO ij = ijpair + 1, o3c%nijpairs
516 17106 : IF (o3c%ijpair(ij)%nklist > 0) THEN
517 14067 : o3c_iterator%ijp_thread(me) = ij
518 14067 : o3c_iterator%k_thread(me) = 1
519 : istat = 0
520 : EXIT
521 : END IF
522 : END DO
523 : ELSE
524 : ! increase klist
525 6839 : o3c_iterator%ijp_thread(me) = ijpair
526 6839 : o3c_iterator%k_thread(me) = klist + 1
527 : istat = 0
528 : END IF
529 :
530 : IF (istat == 0) THEN
531 : ! set last to this thread
532 21018 : o3c_iterator%ijp_last = o3c_iterator%ijp_thread(me)
533 21018 : o3c_iterator%k_last = o3c_iterator%k_thread(me)
534 : ELSE
535 : ! set last to final position
536 112 : o3c_iterator%ijp_last = o3c%nijpairs
537 112 : o3c_iterator%k_last = o3c%ijpair(o3c%nijpairs)%nklist
538 : END IF
539 : !$OMP END CRITICAL(o3c_iterate_critical)
540 :
541 21130 : END FUNCTION o3c_iterate
542 :
543 : ! **************************************************************************************************
544 : !> \brief ...
545 : !> \param o3c_iterator ...
546 : !> \param mepos ...
547 : !> \param integral ...
548 : !> \param tvec ...
549 : !> \param force_i ...
550 : !> \param force_j ...
551 : !> \param force_k ...
552 : ! **************************************************************************************************
553 0 : SUBROUTINE set_o3c_container(o3c_iterator, mepos, integral, tvec, force_i, force_j, force_k)
554 : TYPE(o3c_iterator_type) :: o3c_iterator
555 : INTEGER, OPTIONAL :: mepos
556 : REAL(KIND=dp), DIMENSION(:, :, :), OPTIONAL, &
557 : POINTER :: integral
558 : REAL(KIND=dp), DIMENSION(:, :), OPTIONAL, POINTER :: tvec, force_i, force_j, force_k
559 :
560 : INTEGER :: ij, k, me
561 : TYPE(o3c_container_type), POINTER :: o3c
562 : TYPE(o3c_int_type), POINTER :: ijk
563 : TYPE(o3c_pair_type), POINTER :: ijp
564 :
565 0 : IF (PRESENT(mepos)) THEN
566 0 : me = mepos
567 : ELSE
568 : me = 0
569 : END IF
570 :
571 0 : ij = o3c_iterator%ijp_thread(me)
572 0 : k = o3c_iterator%k_thread(me)
573 :
574 0 : o3c => o3c_iterator%o3c
575 0 : ijp => o3c%ijpair(ij)
576 0 : ijk => ijp%ijk(k)
577 :
578 0 : IF (PRESENT(integral)) ijk%integral => integral
579 0 : IF (PRESENT(tvec)) ijk%tvec => tvec
580 0 : IF (PRESENT(force_i)) ijk%force_i => force_i
581 0 : IF (PRESENT(force_j)) ijk%force_j => force_j
582 0 : IF (PRESENT(force_k)) ijk%force_k => force_k
583 :
584 0 : END SUBROUTINE set_o3c_container
585 :
586 : ! **************************************************************************************************
587 : !> \brief ...
588 : !> \param o3c_vec ...
589 : !> \param nsize ...
590 : ! **************************************************************************************************
591 0 : SUBROUTINE o3c_vec_create(o3c_vec, nsize)
592 : TYPE(o3c_vec_type), DIMENSION(:) :: o3c_vec
593 : INTEGER, DIMENSION(:), INTENT(IN) :: nsize
594 :
595 : INTEGER :: i, m, n
596 :
597 0 : m = SIZE(o3c_vec)
598 0 : CPASSERT(SIZE(nsize) == m)
599 :
600 0 : DO i = 1, m
601 0 : n = nsize(i)
602 0 : ALLOCATE (o3c_vec(i)%v(n))
603 0 : o3c_vec(i)%v = 0.0_dp
604 0 : o3c_vec(i)%n = n
605 : END DO
606 :
607 0 : END SUBROUTINE o3c_vec_create
608 :
609 : ! **************************************************************************************************
610 : !> \brief ...
611 : !> \param o3c_vec ...
612 : ! **************************************************************************************************
613 0 : SUBROUTINE o3c_vec_release(o3c_vec)
614 : TYPE(o3c_vec_type), DIMENSION(:) :: o3c_vec
615 :
616 : INTEGER :: i
617 :
618 0 : DO i = 1, SIZE(o3c_vec)
619 0 : IF (ASSOCIATED(o3c_vec(i)%v)) THEN
620 0 : DEALLOCATE (o3c_vec(i)%v)
621 : END IF
622 : END DO
623 :
624 0 : END SUBROUTINE o3c_vec_release
625 :
626 : ! **************************************************************************************************
627 : !> \brief ...
628 : !> \param o3c_vec ...
629 : !> \param i ...
630 : !> \param vec ...
631 : !> \param n ...
632 : ! **************************************************************************************************
633 0 : SUBROUTINE get_o3c_vec(o3c_vec, i, vec, n)
634 : TYPE(o3c_vec_type), DIMENSION(:) :: o3c_vec
635 : INTEGER, INTENT(IN) :: i
636 : REAL(KIND=dp), DIMENSION(:), OPTIONAL, POINTER :: vec
637 : INTEGER, OPTIONAL :: n
638 :
639 0 : CPASSERT(i > 0 .AND. i <= SIZE(o3c_vec))
640 :
641 0 : IF (PRESENT(vec)) vec => o3c_vec(i)%v
642 0 : IF (PRESENT(n)) n = o3c_vec(i)%n
643 :
644 0 : END SUBROUTINE get_o3c_vec
645 :
646 : ! **************************************************************************************************
647 :
648 0 : END MODULE qs_o3c_types
|
