1!--------------------------------------------------------------------------------------------------! 2! CP2K: A general program to perform molecular dynamics simulations ! 3! Copyright (C) 2000 - 2020 CP2K developers group ! 4!--------------------------------------------------------------------------------------------------! 5 6! ************************************************************************************************** 7!> \brief Define the neighbor list data types and the corresponding functionality 8!> \par History 9!> - cleaned (23.07.2003,MK) 10!> - full refactoring, list iterators (20.10.2010, JGH) 11!> - add get_neighbor_list_set_p, return info for a set of neighborlists 12!> (07.2014,JGH) 13!> \author Matthias Krack (21.06.2000) 14! ************************************************************************************************** 15MODULE qs_neighbor_list_types 16 17 USE kinds, ONLY: dp 18 USE util, ONLY: locate,& 19 sort 20#include "./base/base_uses.f90" 21 22 IMPLICIT NONE 23 24 PRIVATE 25 26! *** Global parameters (in this module) *** 27 28 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_neighbor_list_types' 29 30! *** Definition of the data types for a linked list of neighbors *** 31 32! ************************************************************************************************** 33 TYPE neighbor_node_type 34 PRIVATE 35 TYPE(neighbor_node_type), POINTER :: next_neighbor_node 36 REAL(dp), DIMENSION(3) :: r 37 INTEGER, DIMENSION(3) :: cell 38 INTEGER :: neighbor 39 END TYPE neighbor_node_type 40 41! ************************************************************************************************** 42 TYPE neighbor_list_type 43 PRIVATE 44 TYPE(neighbor_list_type), POINTER :: next_neighbor_list 45 TYPE(neighbor_node_type), POINTER :: first_neighbor_node, & 46 last_neighbor_node 47 INTEGER :: atom, nnode 48 END TYPE neighbor_list_type 49 50! ************************************************************************************************** 51 TYPE neighbor_list_set_type 52 PRIVATE 53 TYPE(neighbor_list_type), POINTER :: first_neighbor_list, & 54 last_neighbor_list 55 INTEGER :: nlist 56 LOGICAL :: symmetric 57 END TYPE neighbor_list_set_type 58 59! ************************************************************************************************** 60 TYPE neighbor_list_p_type 61 TYPE(neighbor_list_type), POINTER :: neighbor_list 62 END TYPE neighbor_list_p_type 63 64! ************************************************************************************************** 65 TYPE neighbor_list_set_p_type 66 TYPE(neighbor_list_set_type), POINTER :: neighbor_list_set 67 INTEGER :: nl_size 68 INTEGER :: nl_start 69 INTEGER :: nl_end 70 TYPE(neighbor_list_task_type), DIMENSION(:), POINTER :: nlist_task 71 END TYPE neighbor_list_set_p_type 72 73! ************************************************************************************************** 74 TYPE list_search_type 75 PRIVATE 76 INTEGER :: nlist 77 INTEGER, DIMENSION(:), POINTER :: atom_list 78 INTEGER, DIMENSION(:), POINTER :: atom_index 79 TYPE(neighbor_list_p_type), & 80 DIMENSION(:), POINTER :: neighbor_list 81 END TYPE list_search_type 82 83! ************************************************************************************************** 84 TYPE neighbor_list_task_type 85 INTEGER :: iatom, jatom, & 86 ikind, jkind, nkind, & 87 ilist, nlist, inode, nnode 88 REAL(KIND=dp), DIMENSION(3) :: r 89 INTEGER, DIMENSION(3) :: cell 90 TYPE(neighbor_list_task_type), & 91 POINTER :: next ! Pointer for forming a linked list of tasks 92 END TYPE neighbor_list_task_type 93 94 INTERFACE nl_sub_iterate 95 MODULE PROCEDURE nl_sub_iterate 96 MODULE PROCEDURE nl_sub_iterate_ref 97 END INTERFACE 98 99! ************************************************************************************************** 100! Neighbor List Iterator 101! ************************************************************************************************** 102 TYPE neighbor_list_iterator_type 103 PRIVATE 104 INTEGER :: ikind, jkind, ilist, inode 105 INTEGER :: nkind, nlist, nnode 106 INTEGER :: iatom, jatom 107 TYPE(neighbor_list_set_p_type), & 108 DIMENSION(:), POINTER :: nl 109 TYPE(neighbor_list_type), POINTER :: neighbor_list 110 TYPE(neighbor_node_type), POINTER :: neighbor_node 111 TYPE(list_search_type), & 112 DIMENSION(:), POINTER :: list_search 113 END TYPE neighbor_list_iterator_type 114 115 TYPE neighbor_list_iterator_p_type 116 PRIVATE 117 TYPE(neighbor_list_iterator_type), POINTER :: neighbor_list_iterator 118 INTEGER :: last 119 END TYPE neighbor_list_iterator_p_type 120! ************************************************************************************************** 121 122! *** Public data types *** 123 124 PUBLIC :: neighbor_list_p_type, & 125 neighbor_list_set_type, & 126 neighbor_list_set_p_type, & 127 neighbor_list_task_type 128 129! *** Public subroutines *** 130 131 PUBLIC :: add_neighbor_list, & 132 add_neighbor_node, & 133 allocate_neighbor_list_set, & 134 deallocate_neighbor_list_set, & 135 release_neighbor_list_sets, & 136 get_iterator_task, & 137 get_neighbor_list_set, & 138 get_neighbor_list_set_p 139 140! *** Iterator functions and types *** 141 142 PUBLIC :: neighbor_list_iterator_p_type, & 143 neighbor_list_iterator_create, & 144 neighbor_list_iterator_release, & 145 neighbor_list_iterate, & 146 nl_set_sub_iterator, & 147 nl_sub_iterate, & 148 get_iterator_info 149 150CONTAINS 151 152! ************************************************************************************************** 153!> \brief Neighbor list iterator functions 154!> \param iterator_set ... 155!> \param nl ... 156!> \param search ... 157!> \param nthread ... 158!> \date 28.07.2010 159!> \author jhu 160!> \version 1.0 161! ************************************************************************************************** 162 SUBROUTINE neighbor_list_iterator_create(iterator_set, nl, search, nthread) 163 TYPE(neighbor_list_iterator_p_type), & 164 DIMENSION(:), POINTER :: iterator_set 165 TYPE(neighbor_list_set_p_type), DIMENSION(:), & 166 POINTER :: nl 167 LOGICAL, INTENT(IN), OPTIONAL :: search 168 INTEGER, INTENT(IN), OPTIONAL :: nthread 169 170 CHARACTER(len=*), PARAMETER :: routineN = 'neighbor_list_iterator_create', & 171 routineP = moduleN//':'//routineN 172 173 INTEGER :: iatom, il, ilist, mthread, nlist 174 TYPE(list_search_type), DIMENSION(:), POINTER :: list_search 175 TYPE(neighbor_list_iterator_type), POINTER :: iterator 176 TYPE(neighbor_list_type), POINTER :: neighbor_list 177 178 mthread = 1 179 IF (PRESENT(nthread)) mthread = nthread 180 181 ALLOCATE (iterator_set(0:mthread - 1)) 182 183 DO il = 0, mthread - 1 184 ALLOCATE (iterator_set(il)%neighbor_list_iterator) 185 186 iterator => iterator_set(il)%neighbor_list_iterator 187 188 iterator%nl => nl 189 190 iterator%ikind = 0 191 iterator%jkind = 0 192 iterator%nkind = NINT(SQRT(REAL(SIZE(nl), dp))) 193 194 iterator%ilist = 0 195 iterator%nlist = 0 196 iterator%inode = 0 197 iterator%nnode = 0 198 199 iterator%iatom = 0 200 iterator%jatom = 0 201 202 NULLIFY (iterator%neighbor_list) 203 NULLIFY (iterator%neighbor_node) 204 NULLIFY (iterator%list_search) 205 END DO 206 207 iterator_set(:)%last = 0 208 209 IF (PRESENT(search)) THEN 210 IF (search) THEN 211 ALLOCATE (list_search(SIZE(nl))) 212 DO il = 1, SIZE(nl) 213 IF (ASSOCIATED(nl(il)%neighbor_list_set)) THEN 214 CALL get_neighbor_list_set(neighbor_list_set=nl(il)%neighbor_list_set, nlist=nlist) 215 list_search(il)%nlist = nlist 216 ALLOCATE (list_search(il)%atom_list(nlist)) 217 ALLOCATE (list_search(il)%atom_index(nlist)) 218 ALLOCATE (list_search(il)%neighbor_list(nlist)) 219 220 NULLIFY (neighbor_list) 221 DO ilist = 1, nlist 222 IF (.NOT. ASSOCIATED(neighbor_list)) THEN 223 neighbor_list => first_list(nl(il)%neighbor_list_set) 224 ELSE 225 neighbor_list => neighbor_list%next_neighbor_list 226 END IF 227 CALL get_neighbor_list(neighbor_list=neighbor_list, atom=iatom) 228 list_search(il)%atom_list(ilist) = iatom 229 list_search(il)%neighbor_list(ilist)%neighbor_list => neighbor_list 230 END DO 231 CALL sort(list_search(il)%atom_list, nlist, list_search(il)%atom_index) 232 233 ELSE 234 list_search(il)%nlist = -1 235 NULLIFY (list_search(il)%atom_list, list_search(il)%atom_index, list_search(il)%neighbor_list) 236 END IF 237 END DO 238 DO il = 0, mthread - 1 239 iterator => iterator_set(il)%neighbor_list_iterator 240 iterator%list_search => list_search 241 END DO 242 END IF 243 END IF 244 245 END SUBROUTINE neighbor_list_iterator_create 246 247! ************************************************************************************************** 248!> \brief ... 249!> \param iterator_set ... 250! ************************************************************************************************** 251 SUBROUTINE neighbor_list_iterator_release(iterator_set) 252 TYPE(neighbor_list_iterator_p_type), & 253 DIMENSION(:), POINTER :: iterator_set 254 255 CHARACTER(len=*), PARAMETER :: routineN = 'neighbor_list_iterator_release', & 256 routineP = moduleN//':'//routineN 257 258 INTEGER :: il, mthread 259 TYPE(neighbor_list_iterator_type), POINTER :: iterator 260 261!all threads have the same search list 262 263 iterator => iterator_set(0)%neighbor_list_iterator 264 IF (ASSOCIATED(iterator%list_search)) THEN 265 DO il = 1, SIZE(iterator%list_search) 266 IF (iterator%list_search(il)%nlist >= 0) THEN 267 DEALLOCATE (iterator%list_search(il)%atom_list) 268 DEALLOCATE (iterator%list_search(il)%atom_index) 269 DEALLOCATE (iterator%list_search(il)%neighbor_list) 270 END IF 271 END DO 272 DEALLOCATE (iterator%list_search) 273 END IF 274 275 mthread = SIZE(iterator_set) 276 DO il = 0, mthread - 1 277 DEALLOCATE (iterator_set(il)%neighbor_list_iterator) 278 END DO 279 DEALLOCATE (iterator_set) 280 281 END SUBROUTINE neighbor_list_iterator_release 282 283! ************************************************************************************************** 284!> \brief ... 285!> \param iterator_set ... 286!> \param ikind ... 287!> \param jkind ... 288!> \param iatom ... 289!> \param mepos ... 290! ************************************************************************************************** 291 SUBROUTINE nl_set_sub_iterator(iterator_set, ikind, jkind, iatom, mepos) 292 TYPE(neighbor_list_iterator_p_type), & 293 DIMENSION(:), POINTER :: iterator_set 294 INTEGER, INTENT(IN) :: ikind, jkind, iatom 295 INTEGER, INTENT(IN), OPTIONAL :: mepos 296 297 CHARACTER(len=*), PARAMETER :: routineN = 'nl_set_sub_iterator', & 298 routineP = moduleN//':'//routineN 299 300 INTEGER :: i, ij, ilist, me, nlist, nnode 301 TYPE(list_search_type), POINTER :: list_search 302 TYPE(neighbor_list_iterator_type), POINTER :: iterator 303 TYPE(neighbor_list_type), POINTER :: neighbor_list 304 305 IF (PRESENT(mepos)) THEN 306 me = mepos 307 ELSE 308 me = 0 309 END IF 310 311 ! Set up my thread-local iterator for the list of iatom / jkind nodes 312 313 iterator => iterator_set(me)%neighbor_list_iterator 314 ij = ikind + iterator%nkind*(jkind - 1) 315 IF (ASSOCIATED(iterator%list_search)) THEN 316 list_search => iterator%list_search(ij) 317 nlist = list_search%nlist 318 ilist = 0 319 NULLIFY (neighbor_list) 320 IF (nlist > 0) THEN 321 i = locate(list_search%atom_list, iatom) 322 i = list_search%atom_index(i) 323 IF (i > 0) neighbor_list => list_search%neighbor_list(i)%neighbor_list 324 ilist = i 325 END IF 326 IF (ASSOCIATED(neighbor_list)) THEN 327 CALL get_neighbor_list(neighbor_list=neighbor_list, nnode=nnode) 328 ELSE 329 nnode = 0 330 END IF 331 ELSE 332 CPABORT("") 333 END IF 334 335 iterator%ikind = ikind 336 iterator%jkind = jkind 337 338 iterator%ilist = ilist 339 iterator%nlist = nlist 340 iterator%inode = 0 341 iterator%nnode = nnode 342 343 iterator%iatom = iatom 344 iterator%jatom = 0 345 346 iterator%neighbor_list => neighbor_list 347 NULLIFY (iterator%neighbor_node) 348 349 END SUBROUTINE nl_set_sub_iterator 350 351! ************************************************************************************************** 352!> \brief ... 353!> \param iterator_set ... 354!> \param mepos ... 355!> \return ... 356! ************************************************************************************************** 357 FUNCTION neighbor_list_iterate(iterator_set, mepos) RESULT(istat) 358 TYPE(neighbor_list_iterator_p_type), & 359 DIMENSION(:), POINTER :: iterator_set 360 INTEGER, OPTIONAL :: mepos 361 INTEGER :: istat 362 363 CHARACTER(len=*), PARAMETER :: routineN = 'neighbor_list_iterate', & 364 routineP = moduleN//':'//routineN 365 366 INTEGER :: iab, last, me 367 TYPE(neighbor_list_iterator_type), POINTER :: iterator 368 TYPE(neighbor_list_set_p_type), DIMENSION(:), & 369 POINTER :: nl 370 371 IF (SIZE(iterator_set) .NE. 1 .AND. .NOT. PRESENT(mepos)) & 372 CPABORT("Parallel iterator calls must include 'mepos'") 373 374 IF (PRESENT(mepos)) THEN 375 me = mepos 376 ELSE 377 me = 0 378 END IF 379 380 istat = 0 381 382!$OMP CRITICAL(neighbour_list_iterate_critical) 383 last = iterator_set(0)%last 384 IF (last /= me) THEN 385 iterator_set(me)%neighbor_list_iterator = iterator_set(last)%neighbor_list_iterator 386 END IF 387 iterator => iterator_set(me)%neighbor_list_iterator 388 nl => iterator%nl 389 390 IF (iterator%inode < iterator%nnode) THEN 391 ! we can be sure that there is another node in this list 392 iterator%inode = iterator%inode + 1 393 iterator%neighbor_node => iterator%neighbor_node%next_neighbor_node 394 ELSE 395 iab = MAX(iterator%ikind + iterator%nkind*(iterator%jkind - 1), 0) 396 kindloop: DO ! look for the next list with nnode /= 0 397 listloop: DO 398 IF (iterator%ilist >= iterator%nlist) EXIT listloop 399 iterator%ilist = iterator%ilist + 1 400 IF (ASSOCIATED(iterator%neighbor_list)) THEN 401 iterator%neighbor_list => iterator%neighbor_list%next_neighbor_list 402 ELSE 403 iterator%neighbor_list => first_list(nl(iab)%neighbor_list_set) 404 END IF 405 CALL get_neighbor_list(neighbor_list=iterator%neighbor_list, atom=iterator%iatom, & 406 nnode=iterator%nnode) 407 IF (iterator%nnode > 0) EXIT kindloop 408 END DO listloop 409 IF (iab >= iterator%nkind**2) THEN 410 istat = 1 411 EXIT kindloop 412 ELSE 413 iab = iab + 1 414 iterator%jkind = (iab - 1)/iterator%nkind + 1 415 iterator%ikind = iab - iterator%nkind*(iterator%jkind - 1) 416 iterator%ilist = 0 417 IF (.NOT. ASSOCIATED(nl(iab)%neighbor_list_set)) THEN 418 iterator%ilist = 0 419 iterator%nlist = 0 420 ELSE 421 CALL get_neighbor_list_set(neighbor_list_set= & 422 nl(iab)%neighbor_list_set, nlist=iterator%nlist) 423 iterator%ilist = 0 424 END IF 425 NULLIFY (iterator%neighbor_list) 426 END IF 427 END DO kindloop 428 IF (istat == 0) THEN 429 iterator%inode = 1 430 iterator%neighbor_node => first_node(iterator%neighbor_list) 431 END IF 432 END IF 433 IF (istat == 0) THEN 434 CALL get_neighbor_node(neighbor_node=iterator%neighbor_node, neighbor=iterator%jatom) 435 END IF 436 437 ! mark the last iterator updated 438 iterator_set(:)%last = me 439!$OMP END CRITICAL(neighbour_list_iterate_critical) 440 441 END FUNCTION neighbor_list_iterate 442 443! ************************************************************************************************** 444!> \brief ... 445!> \param iterator_set ... 446!> \param mepos ... 447!> \return ... 448! ************************************************************************************************** 449 FUNCTION nl_sub_iterate(iterator_set, mepos) RESULT(istat) 450 TYPE(neighbor_list_iterator_p_type), & 451 DIMENSION(:), POINTER :: iterator_set 452 INTEGER, INTENT(IN), OPTIONAL :: mepos 453 INTEGER :: istat 454 455 CHARACTER(len=*), PARAMETER :: routineN = 'nl_sub_iterate', routineP = moduleN//':'//routineN 456 457 INTEGER :: me 458 TYPE(neighbor_list_iterator_type), POINTER :: iterator 459 460 ! Each thread's sub-iterator are independent, no need to synchronise with other threads 461 462 IF (PRESENT(mepos)) THEN 463 me = mepos 464 ELSE 465 me = 0 466 END IF 467 468 istat = 0 469 470 iterator => iterator_set(me)%neighbor_list_iterator 471 472 IF (ASSOCIATED(iterator%neighbor_list)) THEN 473 IF (iterator%inode >= iterator%nnode) THEN 474 ! end of loop 475 istat = 1 476 ELSEIF (iterator%inode == 0) THEN 477 iterator%inode = 1 478 iterator%neighbor_node => first_node(iterator%neighbor_list) 479 ELSEIF (iterator%inode > 0) THEN 480 ! we can be sure that there is another node in this list 481 iterator%inode = iterator%inode + 1 482 iterator%neighbor_node => iterator%neighbor_node%next_neighbor_node 483 ELSE 484 CPABORT("wrong") 485 END IF 486 ELSE 487 ! no list available 488 istat = 1 489 END IF 490 IF (istat == 0) THEN 491 CALL get_neighbor_node(neighbor_node=iterator%neighbor_node, neighbor=iterator%jatom) 492 END IF 493 494 END FUNCTION nl_sub_iterate 495 496! ************************************************************************************************** 497!> \brief wrap nl_sub_iterate s.t. external loop over kinds and calls to nl_set_sub_iterator 498!> are no longer needed. This fixes first atom of iter_sub to second atom of iter_ref. 499!> \param iter_sub ... 500!> \param iter_ref ... 501!> \param mepos ... 502!> \return ... 503! ************************************************************************************************** 504 RECURSIVE FUNCTION nl_sub_iterate_ref(iter_sub, iter_ref, mepos) RESULT(iter_stat) 505 TYPE(neighbor_list_iterator_p_type), & 506 DIMENSION(:), POINTER :: iter_sub, iter_ref 507 INTEGER, INTENT(IN), OPTIONAL :: mepos 508 INTEGER :: iter_stat 509 510 INTEGER :: atom_ref, kind_ref, kind_sub, me, nkind 511 TYPE(neighbor_list_iterator_type), POINTER :: iterator 512 513 IF (PRESENT(mepos)) THEN 514 me = mepos 515 ELSE 516 me = 0 517 END IF 518 519 iterator => iter_sub(me)%neighbor_list_iterator 520 kind_sub = iterator%jkind 521 522 CALL get_iterator_info(iter_ref, jatom=atom_ref, jkind=kind_ref) 523 524 IF (iterator%inode == 0) THEN 525 CALL nl_set_sub_iterator(iter_sub, kind_ref, MAX(kind_sub, 1), atom_ref) 526 ENDIF 527 iter_stat = nl_sub_iterate(iter_sub) 528 IF (iter_stat == 0) RETURN 529 530 nkind = iterator%nkind 531 532 IF (kind_sub == nkind) THEN 533 CALL nl_set_sub_iterator(iter_sub, kind_ref, 1, atom_ref) 534 RETURN 535 ELSE 536 kind_sub = kind_sub + 1 537 CALL nl_set_sub_iterator(iter_sub, kind_ref, kind_sub, atom_ref) 538 iter_stat = nl_sub_iterate_ref(iter_sub, iter_ref) 539 ENDIF 540 541 END FUNCTION 542 543! ************************************************************************************************** 544!> \brief ... 545!> \param iterator_set ... 546!> \param mepos ... 547!> \param ikind ... 548!> \param jkind ... 549!> \param nkind ... 550!> \param ilist ... 551!> \param nlist ... 552!> \param inode ... 553!> \param nnode ... 554!> \param iatom ... 555!> \param jatom ... 556!> \param r ... 557!> \param cell ... 558! ************************************************************************************************** 559 SUBROUTINE get_iterator_info(iterator_set, mepos, & 560 ikind, jkind, nkind, ilist, nlist, inode, nnode, iatom, jatom, r, cell) 561 TYPE(neighbor_list_iterator_p_type), & 562 DIMENSION(:), POINTER :: iterator_set 563 INTEGER, OPTIONAL :: mepos, ikind, jkind, nkind, ilist, & 564 nlist, inode, nnode, iatom, jatom 565 REAL(dp), DIMENSION(3), OPTIONAL :: r 566 INTEGER, DIMENSION(3), OPTIONAL :: cell 567 568 CHARACTER(len=*), PARAMETER :: routineN = 'get_iterator_info', & 569 routineP = moduleN//':'//routineN 570 571 INTEGER :: me 572 TYPE(neighbor_list_iterator_type), POINTER :: iterator 573 574 IF (SIZE(iterator_set) .NE. 1 .AND. .NOT. PRESENT(mepos)) & 575 CPABORT("Parallel iterator calls must include 'mepos'") 576 577 IF (PRESENT(mepos)) THEN 578 me = mepos 579 ELSE 580 me = 0 581 END IF 582 iterator => iterator_set(me)%neighbor_list_iterator 583 584 IF (PRESENT(ikind)) ikind = iterator%ikind 585 IF (PRESENT(jkind)) jkind = iterator%jkind 586 IF (PRESENT(nkind)) nkind = iterator%nkind 587 IF (PRESENT(ilist)) ilist = iterator%ilist 588 IF (PRESENT(nlist)) nlist = iterator%nlist 589 IF (PRESENT(inode)) inode = iterator%inode 590 IF (PRESENT(nnode)) nnode = iterator%nnode 591 IF (PRESENT(iatom)) iatom = iterator%iatom 592 IF (PRESENT(jatom)) jatom = iterator%jatom 593 IF (PRESENT(r)) THEN 594 CALL get_neighbor_node(neighbor_node=iterator%neighbor_node, r=r) 595 END IF 596 IF (PRESENT(cell)) THEN 597 CALL get_neighbor_node(neighbor_node=iterator%neighbor_node, cell=cell) 598 END IF 599 600 END SUBROUTINE get_iterator_info 601 602! ************************************************************************************************** 603!> \brief Captures the current state of the iterator in a neighbor_list_task_type 604!> \param iterator_set the iterator / array of iterators (for multiple threads) 605!> \param task the task structure which is returned 606!> \param mepos OpenMP thread index 607! ************************************************************************************************** 608 SUBROUTINE get_iterator_task(iterator_set, task, mepos) 609 TYPE(neighbor_list_iterator_p_type), & 610 DIMENSION(:), POINTER :: iterator_set 611 TYPE(neighbor_list_task_type), INTENT(OUT) :: task 612 INTEGER, OPTIONAL :: mepos 613 614 IF (PRESENT(mepos)) THEN 615 CALL get_iterator_info(iterator_set, mepos=mepos, ikind=task%ikind, jkind=task%jkind, & 616 nkind=task%nkind, & 617 ilist=task%ilist, nlist=task%nlist, & 618 inode=task%inode, nnode=task%nnode, & 619 iatom=task%iatom, jatom=task%jatom, & 620 r=task%r, cell=task%cell) 621 ELSE 622 CALL get_iterator_info(iterator_set, ikind=task%ikind, jkind=task%jkind, & 623 nkind=task%nkind, & 624 ilist=task%ilist, nlist=task%nlist, & 625 inode=task%inode, nnode=task%nnode, & 626 iatom=task%iatom, jatom=task%jatom, & 627 r=task%r, cell=task%cell) 628 ENDIF 629 630 NULLIFY (task%next) 631 632 END SUBROUTINE get_iterator_task 633 634! ************************************************************************************************** 635!> \brief Add a new neighbor list to a neighbor list set. 636!> \param neighbor_list_set ... 637!> \param atom ... 638!> \param neighbor_list ... 639!> \date 13.09.2000 640!> \author MK 641!> \version 1.0 642! ************************************************************************************************** 643 SUBROUTINE add_neighbor_list(neighbor_list_set, atom, neighbor_list) 644 645 TYPE(neighbor_list_set_type), POINTER :: neighbor_list_set 646 INTEGER, INTENT(IN) :: atom 647 TYPE(neighbor_list_type), POINTER :: neighbor_list 648 649 CHARACTER(len=*), PARAMETER :: routineN = 'add_neighbor_list', & 650 routineP = moduleN//':'//routineN 651 652 TYPE(neighbor_list_type), POINTER :: new_neighbor_list 653 654 IF (ASSOCIATED(neighbor_list_set)) THEN 655 656 IF (ASSOCIATED(neighbor_list_set%last_neighbor_list)) THEN 657 658 new_neighbor_list => & 659 neighbor_list_set%last_neighbor_list%next_neighbor_list 660 661 IF (.NOT. ASSOCIATED(new_neighbor_list)) THEN 662 663! *** Allocate a new neighbor list *** 664 665 ALLOCATE (new_neighbor_list) 666 667 NULLIFY (new_neighbor_list%next_neighbor_list) 668 NULLIFY (new_neighbor_list%first_neighbor_node) 669 670! *** Link the new neighbor list to the neighbor list set *** 671 672 neighbor_list_set%last_neighbor_list%next_neighbor_list => new_neighbor_list 673 674 END IF 675 676 ELSE 677 678 new_neighbor_list => neighbor_list_set%first_neighbor_list 679 680 IF (.NOT. ASSOCIATED(new_neighbor_list)) THEN 681 682! *** Allocate a new first neighbor list *** 683 684 ALLOCATE (new_neighbor_list) 685 686 NULLIFY (new_neighbor_list%next_neighbor_list) 687 NULLIFY (new_neighbor_list%first_neighbor_node) 688 689! *** Link the new first neighbor list to the neighbor list set *** 690 691 neighbor_list_set%first_neighbor_list => new_neighbor_list 692 693 END IF 694 695 END IF 696 697! *** Store the data set of the new neighbor list *** 698 699 NULLIFY (new_neighbor_list%last_neighbor_node) 700 new_neighbor_list%atom = atom 701 new_neighbor_list%nnode = 0 702 703! *** Update the pointer to the last neighbor *** 704! *** list of the neighbor list set *** 705 706 neighbor_list_set%last_neighbor_list => new_neighbor_list 707 708! *** Increment the neighbor list counter *** 709 710 neighbor_list_set%nlist = neighbor_list_set%nlist + 1 711 712! *** Return a pointer to the new neighbor list *** 713 714 neighbor_list => new_neighbor_list 715 716 ELSE 717 718 CPABORT("The requested neighbor list set is not associated") 719 720 END IF 721 722 END SUBROUTINE add_neighbor_list 723 724! ************************************************************************************************** 725!> \brief Add a new neighbor list node to a neighbor list. 726!> \param neighbor_list ... 727!> \param neighbor ... 728!> \param cell ... 729!> \param r ... 730!> \param exclusion_list ... 731!> \param nkind ... 732!> \date 23.06.2000 733!> \author MK 734!> \version 1.0 735! ************************************************************************************************** 736 SUBROUTINE add_neighbor_node(neighbor_list, neighbor, cell, r, exclusion_list, nkind) 737 738 TYPE(neighbor_list_type), POINTER :: neighbor_list 739 INTEGER, INTENT(IN) :: neighbor 740 INTEGER, DIMENSION(3), INTENT(IN) :: cell 741 REAL(dp), DIMENSION(3), INTENT(IN) :: r 742 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: exclusion_list 743 INTEGER, INTENT(IN), OPTIONAL :: nkind 744 745 CHARACTER(len=*), PARAMETER :: routineN = 'add_neighbor_node', & 746 routineP = moduleN//':'//routineN 747 748 INTEGER :: iatom, my_nkind 749 TYPE(neighbor_node_type), POINTER :: new_neighbor_node 750 751 IF (ASSOCIATED(neighbor_list)) THEN 752 753! *** Check for exclusions *** 754 755 IF (PRESENT(exclusion_list)) THEN 756 IF (ASSOCIATED(exclusion_list)) THEN 757 DO iatom = 1, SIZE(exclusion_list) 758 IF (exclusion_list(iatom) == 0) EXIT 759 IF (exclusion_list(iatom) == neighbor) RETURN 760 END DO 761 END IF 762 END IF 763 764 my_nkind = 0 765 IF (PRESENT(nkind)) my_nkind = nkind 766 767 IF (ASSOCIATED(neighbor_list%last_neighbor_node)) THEN 768 769 new_neighbor_node => neighbor_list%last_neighbor_node%next_neighbor_node 770 771 IF (.NOT. ASSOCIATED(new_neighbor_node)) THEN 772 773! *** Allocate a new neighbor node *** 774 775 ALLOCATE (new_neighbor_node) 776 777 NULLIFY (new_neighbor_node%next_neighbor_node) 778 779! *** Link the new neighbor node to the neighbor list *** 780 781 neighbor_list%last_neighbor_node%next_neighbor_node => new_neighbor_node 782 783 END IF 784 785 ELSE 786 787 new_neighbor_node => neighbor_list%first_neighbor_node 788 789 IF (.NOT. ASSOCIATED(new_neighbor_node)) THEN 790 791! *** Allocate a new first neighbor node *** 792 793 ALLOCATE (new_neighbor_node) 794 795 NULLIFY (new_neighbor_node%next_neighbor_node) 796 797! *** Link the new first neighbor node to the neighbor list *** 798 799 neighbor_list%first_neighbor_node => new_neighbor_node 800 801 END IF 802 803 END IF 804 805! *** Store the data set of the new neighbor *** 806 807 new_neighbor_node%neighbor = neighbor 808 new_neighbor_node%cell(:) = cell(:) 809 new_neighbor_node%r(:) = r(:) 810 811! *** Update the pointer to the last neighbor node of the neighbor list *** 812 813 neighbor_list%last_neighbor_node => new_neighbor_node 814 815! *** Increment the neighbor node counter *** 816 817 neighbor_list%nnode = neighbor_list%nnode + 1 818 819 ELSE 820 821 CPABORT("The requested neighbor list is not associated") 822 823 END IF 824 825 END SUBROUTINE add_neighbor_node 826 827! ************************************************************************************************** 828!> \brief Allocate and initialize a set of neighbor lists. 829!> \param neighbor_list_set ... 830!> \param symmetric ... 831!> \date 23.06.2000 832!> \author MK 833!> \version 1.0 834! ************************************************************************************************** 835 SUBROUTINE allocate_neighbor_list_set(neighbor_list_set, symmetric) 836 837 TYPE(neighbor_list_set_type), POINTER :: neighbor_list_set 838 LOGICAL, INTENT(IN) :: symmetric 839 840 CHARACTER(len=*), PARAMETER :: routineN = 'allocate_neighbor_list_set', & 841 routineP = moduleN//':'//routineN 842 843! *** Deallocate the old neighbor list set *** 844 845 IF (ASSOCIATED(neighbor_list_set)) THEN 846 CALL deallocate_neighbor_list_set(neighbor_list_set) 847 END IF 848 849! *** Allocate a set of neighbor lists *** 850 851 ALLOCATE (neighbor_list_set) 852 853 NULLIFY (neighbor_list_set%first_neighbor_list) 854 855! *** Initialize the pointers to the first neighbor list *** 856 857 CALL init_neighbor_list_set(neighbor_list_set, symmetric) 858 859 END SUBROUTINE allocate_neighbor_list_set 860 861! ************************************************************************************************** 862!> \brief Deallocate a neighbor list. 863!> \param neighbor_list ... 864!> \date 20.09.2002 865!> \author MK 866!> \version 1.0 867! ************************************************************************************************** 868 SUBROUTINE deallocate_neighbor_list(neighbor_list) 869 870 TYPE(neighbor_list_type), POINTER :: neighbor_list 871 872 CHARACTER(len=*), PARAMETER :: routineN = 'deallocate_neighbor_list', & 873 routineP = moduleN//':'//routineN 874 875 TYPE(neighbor_node_type), POINTER :: neighbor_node, next_neighbor_node 876 877 IF (ASSOCIATED(neighbor_list)) THEN 878 879 neighbor_node => neighbor_list%first_neighbor_node 880 881 DO WHILE (ASSOCIATED(neighbor_node)) 882 next_neighbor_node => neighbor_node%next_neighbor_node 883 DEALLOCATE (neighbor_node) 884 neighbor_node => next_neighbor_node 885 END DO 886 887 DEALLOCATE (neighbor_list) 888 889 END IF 890 891 END SUBROUTINE deallocate_neighbor_list 892 893! ************************************************************************************************** 894!> \brief Deallocate a neighbor list set. 895!> \param neighbor_list_set ... 896!> \date 03.11.2000 897!> \author MK 898!> \version 1.0 899! ************************************************************************************************** 900 SUBROUTINE deallocate_neighbor_list_set(neighbor_list_set) 901 TYPE(neighbor_list_set_type), POINTER :: neighbor_list_set 902 903 CHARACTER(len=*), PARAMETER :: routineN = 'deallocate_neighbor_list_set', & 904 routineP = moduleN//':'//routineN 905 906 TYPE(neighbor_list_type), POINTER :: neighbor_list, next_neighbor_list 907 908 IF (ASSOCIATED(neighbor_list_set)) THEN 909 910 neighbor_list => neighbor_list_set%first_neighbor_list 911 912 DO WHILE (ASSOCIATED(neighbor_list)) 913 next_neighbor_list => neighbor_list%next_neighbor_list 914 CALL deallocate_neighbor_list(neighbor_list) 915 neighbor_list => next_neighbor_list 916 END DO 917 918 DEALLOCATE (neighbor_list_set) 919 920 END IF 921 922 END SUBROUTINE deallocate_neighbor_list_set 923 924! ************************************************************************************************** 925!> \brief Return a pointer to the first neighbor list of a neighbor list set. 926!> \param neighbor_list_set ... 927!> \return ... 928!> \date 13.09.2000 929!> \author MK 930!> \version 1.0 931! ************************************************************************************************** 932 FUNCTION first_list(neighbor_list_set) RESULT(first_neighbor_list) 933 934 TYPE(neighbor_list_set_type), POINTER :: neighbor_list_set 935 TYPE(neighbor_list_type), POINTER :: first_neighbor_list 936 937 first_neighbor_list => neighbor_list_set%first_neighbor_list 938 939 END FUNCTION first_list 940 941! ************************************************************************************************** 942!> \brief Return a pointer to the first neighbor node of a neighbor list. 943!> \param neighbor_list ... 944!> \return ... 945!> \date 23.06.2000, 946!> \author MK 947!> \version 1.0 948! ************************************************************************************************** 949 FUNCTION first_node(neighbor_list) RESULT(first_neighbor_node) 950 951 TYPE(neighbor_list_type), POINTER :: neighbor_list 952 TYPE(neighbor_node_type), POINTER :: first_neighbor_node 953 954 first_neighbor_node => neighbor_list%first_neighbor_node 955 956 END FUNCTION first_node 957 958! ************************************************************************************************** 959!> \brief Return the requested data of a neighbor list. 960!> \param neighbor_list ... 961!> \param atom ... 962!> \param nnode ... 963!> \date 13.09.2000 964!> \author MK 965!> \version 1.0 966! ************************************************************************************************** 967 SUBROUTINE get_neighbor_list(neighbor_list, atom, nnode) 968 969 TYPE(neighbor_list_type), POINTER :: neighbor_list 970 INTEGER, INTENT(OUT), OPTIONAL :: atom, nnode 971 972 CHARACTER(len=*), PARAMETER :: routineN = 'get_neighbor_list', & 973 routineP = moduleN//':'//routineN 974 975 IF (ASSOCIATED(neighbor_list)) THEN 976 977 IF (PRESENT(atom)) atom = neighbor_list%atom 978 IF (PRESENT(nnode)) nnode = neighbor_list%nnode 979 980 ELSE 981 982 CPABORT("The requested neighbor list is not associated") 983 984 END IF 985 986 END SUBROUTINE get_neighbor_list 987 988! ************************************************************************************************** 989!> \brief Return the components of a neighbor list set. 990!> \param neighbor_list_set ... 991!> \param nlist ... 992!> \param symmetric ... 993!> \date 10.11.2000 994!> \author MK 995!> \version 1.0 996! ************************************************************************************************** 997 SUBROUTINE get_neighbor_list_set(neighbor_list_set, nlist, symmetric) 998 999 TYPE(neighbor_list_set_type), POINTER :: neighbor_list_set 1000 INTEGER, INTENT(OUT), OPTIONAL :: nlist 1001 LOGICAL, INTENT(OUT), OPTIONAL :: symmetric 1002 1003 CHARACTER(len=*), PARAMETER :: routineN = 'get_neighbor_list_set', & 1004 routineP = moduleN//':'//routineN 1005 1006 IF (ASSOCIATED(neighbor_list_set)) THEN 1007 1008 IF (PRESENT(nlist)) nlist = neighbor_list_set%nlist 1009 IF (PRESENT(symmetric)) symmetric = neighbor_list_set%symmetric 1010 1011 ELSE 1012 1013 CPABORT("The requested neighbor list set is not associated") 1014 1015 END IF 1016 1017 END SUBROUTINE get_neighbor_list_set 1018 1019! ************************************************************************************************** 1020!> \brief Return the components of the first neighbor list set. 1021!> \param neighbor_list_sets ... 1022!> \param nlist ... 1023!> \param symmetric ... 1024!> \date 07.2014 1025!> \author JGH 1026!> \version 1.0 1027! ************************************************************************************************** 1028 SUBROUTINE get_neighbor_list_set_p(neighbor_list_sets, nlist, symmetric) 1029 1030 TYPE(neighbor_list_set_p_type), DIMENSION(:), & 1031 POINTER :: neighbor_list_sets 1032 INTEGER, INTENT(OUT), OPTIONAL :: nlist 1033 LOGICAL, INTENT(OUT), OPTIONAL :: symmetric 1034 1035 CHARACTER(len=*), PARAMETER :: routineN = 'get_neighbor_list_set_p', & 1036 routineP = moduleN//':'//routineN 1037 1038 INTEGER :: i 1039 TYPE(neighbor_list_set_type), POINTER :: neighbor_list_set 1040 1041 IF (ASSOCIATED(neighbor_list_sets)) THEN 1042 1043 NULLIFY (neighbor_list_set) 1044 DO i = 1, SIZE(neighbor_list_sets) 1045 neighbor_list_set => neighbor_list_sets(i)%neighbor_list_set 1046 IF (ASSOCIATED(neighbor_list_set)) EXIT 1047 END DO 1048 1049 IF (ASSOCIATED(neighbor_list_set)) THEN 1050 IF (PRESENT(nlist)) nlist = neighbor_list_set%nlist 1051 IF (PRESENT(symmetric)) symmetric = neighbor_list_set%symmetric 1052 ELSE 1053 CALL cp_abort(__LOCATION__, "No neighbor list set is associated. "// & 1054 "Did you specify *all* required basis-sets, eg. for ADMM?") 1055 END IF 1056 1057 ELSE 1058 1059 CPABORT("The requested neighbor list sets are not associated") 1060 1061 END IF 1062 1063 END SUBROUTINE get_neighbor_list_set_p 1064 1065! ************************************************************************************************** 1066!> \brief Return the requested data of a neighbor node. 1067!> \param neighbor_node ... 1068!> \param neighbor ... 1069!> \param cell ... 1070!> \param r ... 1071!> \date 23.06.2000 1072!> \author MK 1073!> \version 1.0 1074! ************************************************************************************************** 1075 SUBROUTINE get_neighbor_node(neighbor_node, neighbor, cell, r) 1076 1077 TYPE(neighbor_node_type), POINTER :: neighbor_node 1078 INTEGER, INTENT(OUT), OPTIONAL :: neighbor 1079 INTEGER, DIMENSION(3), INTENT(OUT), OPTIONAL :: cell 1080 REAL(dp), DIMENSION(3), INTENT(OUT), OPTIONAL :: r 1081 1082 CHARACTER(len=*), PARAMETER :: routineN = 'get_neighbor_node', & 1083 routineP = moduleN//':'//routineN 1084 1085 IF (ASSOCIATED(neighbor_node)) THEN 1086 1087 IF (PRESENT(neighbor)) neighbor = neighbor_node%neighbor 1088 IF (PRESENT(r)) r(:) = neighbor_node%r(:) 1089 IF (PRESENT(cell)) cell(:) = neighbor_node%cell(:) 1090 1091 ELSE 1092 1093 CPABORT("The requested neighbor node is not associated") 1094 1095 END IF 1096 1097 END SUBROUTINE get_neighbor_node 1098 1099! ************************************************************************************************** 1100!> \brief Initialize a neighbor list set. Nothing is (de)allocated here. 1101!> This routine is also used to prepare a neighbor list set for 1102!> overwriting. 1103!> \param neighbor_list_set ... 1104!> \param symmetric ... 1105!> \date 20.09.2002 1106!> \author MK 1107!> \version 1.0 1108! ************************************************************************************************** 1109 SUBROUTINE init_neighbor_list_set(neighbor_list_set, symmetric) 1110 1111 TYPE(neighbor_list_set_type), POINTER :: neighbor_list_set 1112 LOGICAL, INTENT(IN) :: symmetric 1113 1114 CHARACTER(len=*), PARAMETER :: routineN = 'init_neighbor_list_set', & 1115 routineP = moduleN//':'//routineN 1116 1117 IF (ASSOCIATED(neighbor_list_set)) THEN 1118 1119 ! *** Initialize the pointers to the last neighbor list *** 1120 NULLIFY (neighbor_list_set%last_neighbor_list) 1121 1122 ! *** Initialize the neighbor list counter *** 1123 neighbor_list_set%nlist = 0 1124 1125 ! *** Initialize the neighbor list build properties 1126 neighbor_list_set%symmetric = symmetric 1127 1128 ELSE 1129 1130 CPABORT("The requested neighbor list set is not associated") 1131 1132 END IF 1133 1134 END SUBROUTINE init_neighbor_list_set 1135 1136! ************************************************************************************************** 1137!> \brief releases an array of neighbor_list_sets 1138!> \param nlists ... 1139!> \author Ole Schuett 1140! ************************************************************************************************** 1141 SUBROUTINE release_neighbor_list_sets(nlists) 1142 TYPE(neighbor_list_set_p_type), DIMENSION(:), & 1143 POINTER :: nlists 1144 1145 INTEGER :: i 1146 1147 IF (ASSOCIATED(nlists)) THEN 1148 DO i = 1, SIZE(nlists) 1149 CALL deallocate_neighbor_list_set(nlists(i)%neighbor_list_set) 1150 END DO 1151 IF (ASSOCIATED(nlists(1)%nlist_task)) THEN 1152 DEALLOCATE (nlists(1)%nlist_task) 1153 END IF 1154 DEALLOCATE (nlists) 1155 END IF 1156 END SUBROUTINE release_neighbor_list_sets 1157 1158END MODULE qs_neighbor_list_types 1159