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