290 LOGICAL,
OPTIONAL :: molecular
293 CHARACTER(len=*),
PARAMETER :: routinen =
'build_qs_neighbor_lists'
295 CHARACTER(LEN=2) :: element_symbol, element_symbol2
296 CHARACTER(LEN=default_string_length) :: print_key_path
297 INTEGER :: handle, hfx_pot, ikind, ingp, iw, jkind, &
298 maxatom, ngp, nkind, zat
299 LOGICAL :: all_potential_present, almo, cneo_potential_present, dftb, do_hfx, dokp, &
300 gth_potential_present, lri_optbas, lrigpw, mic, molecule_only, nddo, paw_atom, &
301 paw_atom_present, rigpw, sgp_potential_present, xtb
302 LOGICAL,
ALLOCATABLE,
DIMENSION(:) :: all_present, aux_fit_present, aux_present, &
303 cneo_present, core_present, default_present, nonbond1_atom, nonbond2_atom, oce_present, &
304 orb_present, ppl_present, ppnl_present, ri_present, xb1_atom, xb2_atom
305 REAL(
dp) :: almo_rcov, almo_rvdw, eps_schwarz, &
306 omega, pdist, rcut, roperator, subcells
307 REAL(
dp),
ALLOCATABLE,
DIMENSION(:) :: all_pot_rad, aux_fit_radius, c_radius, calpha, &
308 core_radius, nuc_orb_radius, oce_radius, orb_radius, ppl_radius, ppnl_radius, ri_radius, &
310 REAL(
dp),
ALLOCATABLE,
DIMENSION(:, :) :: pair_radius, pair_radius_lb
322 nuc_basis_set, orb_basis_set, &
328 sab_cn, sab_cneo, sab_core, sab_gcp, sab_kp, sab_kp_nosym, sab_lrc, sab_orb, sab_scp, &
329 sab_se, sab_tbe, sab_vdw, sab_xb, sab_xtb_nonbond, sab_xtb_pp, sab_xtbe, sac_ae, sac_lri, &
330 sac_ppl, sap_oce, sap_ppnl, soa_list, soo_list
336 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
342 CALL timeset(routinen, handle)
346 NULLIFY (atomic_kind_set, qs_kind_set, cell, neighbor_list_section, &
347 distribution_1d, distribution_2d, gth_potential, sgp_potential, orb_basis_set, &
348 particle_set, molecule_set, dft_control, ks_env)
363 NULLIFY (sab_xtb_nonbond)
371 NULLIFY (sab_kp_nosym)
376 atomic_kind_set=atomic_kind_set, &
377 qs_kind_set=qs_kind_set, &
380 distribution_2d=distribution_2d, &
381 local_particles=distribution_1d, &
382 particle_set=particle_set, &
383 molecule_set=molecule_set, &
384 dft_control=dft_control)
390 last_qs_neighbor_list_id_nr = last_qs_neighbor_list_id_nr + 1
391 CALL set_ks_env(ks_env=ks_env, neighbor_list_id=last_qs_neighbor_list_id_nr)
407 sab_xtb_pp=sab_xtb_pp, &
408 sab_xtb_nonbond=sab_xtb_nonbond, &
413 sab_kp_nosym=sab_kp_nosym, &
416 dokp = (kpoints%nkp > 0)
417 nddo = dft_control%qs_control%semi_empirical
418 dftb = dft_control%qs_control%dftb
419 xtb = dft_control%qs_control%xtb
420 almo = dft_control%qs_control%do_almo_scf
423 lri_optbas = dft_control%qs_control%lri_optbas
426 molecule_only = .false.
427 IF (
PRESENT(molecular)) molecule_only = molecular
437 pdist = dft_control%qs_control%pairlist_radius
444 gth_potential_present=gth_potential_present, &
445 sgp_potential_present=sgp_potential_present, &
446 all_potential_present=all_potential_present, &
447 cneo_potential_present=cneo_potential_present)
452 nkind =
SIZE(atomic_kind_set)
453 ALLOCATE (orb_present(nkind), aux_fit_present(nkind), aux_present(nkind), &
454 default_present(nkind), core_present(nkind))
455 ALLOCATE (orb_radius(nkind), aux_fit_radius(nkind), c_radius(nkind), &
456 core_radius(nkind), calpha(nkind), zeff(nkind))
457 orb_radius(:) = 0.0_dp
458 aux_fit_radius(:) = 0.0_dp
460 core_radius(:) = 0.0_dp
464 ALLOCATE (pair_radius(nkind, nkind))
465 IF (gth_potential_present .OR. sgp_potential_present)
THEN
466 ALLOCATE (ppl_present(nkind), ppl_radius(nkind))
468 ALLOCATE (ppnl_present(nkind), ppnl_radius(nkind))
471 IF (paw_atom_present)
THEN
472 ALLOCATE (oce_present(nkind), oce_radius(nkind))
475 IF (all_potential_present .OR. sgp_potential_present)
THEN
476 ALLOCATE (all_present(nkind), all_pot_rad(nkind))
479 IF (cneo_potential_present)
THEN
480 ALLOCATE (cneo_present(nkind), nuc_orb_radius(nkind))
481 nuc_orb_radius = 0.0_dp
485 ALLOCATE (atom2d(nkind))
486 CALL atom2d_build(atom2d, distribution_1d, distribution_2d, atomic_kind_set, &
487 molecule_set, molecule_only, particle_set=particle_set)
491 CALL get_atomic_kind(atomic_kind_set(ikind), atom_list=atom2d(ikind)%list)
493 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set, basis_type=
"ORB")
495 CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_fit_basis_set, basis_type=
"AUX_FIT")
496 CALL get_qs_kind(qs_kind_set(ikind), basis_set=nuc_basis_set, basis_type=
"NUC")
499 paw_proj_set=paw_proj, &
501 all_potential=all_potential, &
502 gth_potential=gth_potential, &
503 sgp_potential=sgp_potential, &
504 cneo_potential=cneo_potential)
509 CALL get_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom)
511 cutoff=orb_radius(ikind), &
512 defined=orb_present(ikind))
514 IF (
ASSOCIATED(orb_basis_set))
THEN
515 orb_present(ikind) = .true.
516 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, kind_radius=orb_radius(ikind))
518 orb_present(ikind) = .false.
523 aux_present(ikind) = .true.
525 aux_present(ikind) = .false.
528 IF (
ASSOCIATED(aux_fit_basis_set))
THEN
529 aux_fit_present(ikind) = .true.
530 CALL get_gto_basis_set(gto_basis_set=aux_fit_basis_set, kind_radius=aux_fit_radius(ikind))
532 aux_fit_present(ikind) = .false.
535 core_present(ikind) = .false.
536 IF (
ASSOCIATED(cneo_potential) .AND.
ASSOCIATED(nuc_basis_set))
THEN
537 cneo_present(ikind) = .true.
538 CALL get_gto_basis_set(gto_basis_set=nuc_basis_set, kind_radius=nuc_orb_radius(ikind))
540 IF (cneo_potential_present) cneo_present(ikind) = .false.
543 alpha_core_charge=calpha(ikind), &
544 core_charge_radius=core_radius(ikind), &
546 IF (zeff(ikind) /= 0._dp .AND. calpha(ikind) /= 0._dp)
THEN
547 core_present(ikind) = .true.
549 core_present(ikind) = .false.
554 IF (gth_potential_present .OR. sgp_potential_present)
THEN
555 IF (
ASSOCIATED(gth_potential))
THEN
557 ppl_present=ppl_present(ikind), &
558 ppl_radius=ppl_radius(ikind), &
559 ppnl_present=ppnl_present(ikind), &
560 ppnl_radius=ppnl_radius(ikind))
561 ELSE IF (
ASSOCIATED(sgp_potential))
THEN
563 ppl_present=ppl_present(ikind), &
564 ppl_radius=ppl_radius(ikind), &
565 ppnl_present=ppnl_present(ikind), &
566 ppnl_radius=ppnl_radius(ikind))
568 ppl_present(ikind) = .false.
569 ppnl_present(ikind) = .false.
574 IF (paw_atom_present)
THEN
576 oce_present(ikind) = .true.
579 oce_present(ikind) = .false.
584 IF (all_potential_present .OR. sgp_potential_present)
THEN
585 all_present(ikind) = .false.
586 all_pot_rad(ikind) = 0.0_dp
587 IF (
ASSOCIATED(all_potential))
THEN
588 all_present(ikind) = .true.
589 CALL get_potential(potential=all_potential, core_charge_radius=all_pot_rad(ikind))
590 ELSE IF (
ASSOCIATED(sgp_potential))
THEN
591 IF (sgp_potential%ecp_local)
THEN
592 all_present(ikind) = .true.
593 CALL get_potential(potential=sgp_potential, core_charge_radius=all_pot_rad(ikind))
601 IF (pdist < 0.0_dp)
THEN
606 CALL pair_radius_setup(orb_present, orb_present, orb_radius, orb_radius, pair_radius, pdist)
608 mic=mic, subcells=subcells, molecular=molecule_only, nlname=
"sab_orb")
609 CALL set_ks_env(ks_env=ks_env, sab_orb=sab_orb)
611 "/SAB_ORB",
"sab_orb",
"ORBITAL ORBITAL")
616 IF (.NOT. (nddo .OR. dftb .OR. xtb))
THEN
618 mic=mic, symmetric=.false., subcells=subcells, molecular=molecule_only, nlname=
"sab_all")
619 CALL set_ks_env(ks_env=ks_env, sab_all=sab_all)
623 IF (.NOT. (nddo .OR. dftb .OR. xtb))
THEN
624 CALL pair_radius_setup(core_present, core_present, core_radius, core_radius, pair_radius)
625 CALL build_neighbor_lists(sab_core, particle_set, atom2d, cell, pair_radius, subcells=subcells, &
626 operator_type=
"PP", nlname=
"sab_core")
627 CALL set_ks_env(ks_env=ks_env, sab_core=sab_core)
629 "/SAB_CORE",
"sab_core",
"CORE CORE")
642 SELECT CASE (hfx_pot)
654 cpabort(
"HFX potential not available for K-points (NYI)")
657 IF (dft_control%do_admm)
THEN
658 CALL pair_radius_setup(aux_fit_present, aux_fit_present, aux_fit_radius, aux_fit_radius, &
662 ALLOCATE (pair_radius_lb(nkind, nkind))
663 CALL pair_radius_setup(orb_present, orb_present, orb_radius, orb_radius, pair_radius_lb)
666 IF (pair_radius(ikind, jkind) +
cutoff_screen_factor*roperator .LE. pair_radius_lb(ikind, jkind)) &
667 pair_radius(ikind, jkind) = pair_radius_lb(ikind, jkind) - roperator
671 CALL pair_radius_setup(orb_present, orb_present, orb_radius, orb_radius, pair_radius)
675 CALL pair_radius_setup(orb_present, orb_present, orb_radius, orb_radius, pair_radius)
678 subcells=subcells, nlname=
"sab_kp")
683 subcells=subcells, nlname=
"sab_kp_nosym", symmetric=.false.)
684 CALL set_ks_env(ks_env=ks_env, sab_kp_nosym=sab_kp_nosym)
689 IF (gth_potential_present .OR. sgp_potential_present)
THEN
690 IF (any(ppl_present))
THEN
691 CALL pair_radius_setup(orb_present, ppl_present, orb_radius, ppl_radius, pair_radius)
693 subcells=subcells, operator_type=
"ABC", nlname=
"sac_ppl")
694 CALL set_ks_env(ks_env=ks_env, sac_ppl=sac_ppl)
696 "/SAC_PPL",
"sac_ppl",
"ORBITAL GTH-PPL")
698 IF (qs_env%lri_env%ppl_ri)
THEN
700 subcells=subcells, symmetric=.false., operator_type=
"PP", nlname=
"sac_lri")
701 CALL set_ks_env(ks_env=ks_env, sac_lri=sac_lri)
706 IF (any(ppnl_present))
THEN
707 CALL pair_radius_setup(orb_present, ppnl_present, orb_radius, ppnl_radius, pair_radius)
709 subcells=subcells, operator_type=
"ABBA", nlname=
"sap_ppnl")
710 CALL set_ks_env(ks_env=ks_env, sap_ppnl=sap_ppnl)
712 "/SAP_PPNL",
"sap_ppnl",
"ORBITAL GTH-PPNL")
716 IF (paw_atom_present)
THEN
718 IF (any(oce_present))
THEN
719 CALL pair_radius_setup(orb_present, oce_present, orb_radius, oce_radius, pair_radius)
721 subcells=subcells, operator_type=
"ABBA", nlname=
"sap_oce")
722 CALL set_ks_env(ks_env=ks_env, sap_oce=sap_oce)
724 "/SAP_OCE",
"sap_oce",
"ORBITAL(A) PAW-PRJ")
729 IF (.NOT. (nddo .OR. dftb .OR. xtb))
THEN
730 IF (all_potential_present .OR. sgp_potential_present)
THEN
731 CALL pair_radius_setup(orb_present, all_present, orb_radius, all_pot_rad, pair_radius)
733 subcells=subcells, operator_type=
"ABC", nlname=
"sac_ae")
736 "/SAC_AE",
"sac_ae",
"ORBITAL ERFC POTENTIAL")
741 IF (cneo_potential_present)
THEN
742 CALL pair_radius_setup(cneo_present, core_present, nuc_orb_radius, core_radius, pair_radius)
744 subcells=subcells, symmetric=.false., operator_type=
"PP", nlname=
"sab_cneo")
745 CALL set_ks_env(ks_env=ks_env, sab_cneo=sab_cneo)
747 "/SAB_CNEO",
"sab_cneo",
"NUCLEAR ORBITAL ERFC POTENTIAL")
752 default_present = .true.
753 c_radius = dft_control%qs_control%se_control%cutoff_cou
755 CALL pair_radius_setup(default_present, default_present, c_radius, c_radius, pair_radius)
756 IF (dft_control%qs_control%se_control%do_ewald_gks)
THEN
759 subcells=subcells, nlname=
"sab_se")
762 subcells=subcells, nlname=
"sab_se")
766 "/SAB_SE",
"sab_se",
"HARTREE INTERACTIONS")
769 IF ((dft_control%qs_control%se_control%do_ewald) .AND. &
770 (dft_control%qs_control%se_control%integral_screening /=
do_se_is_slater))
THEN
771 c_radius = dft_control%qs_control%se_control%cutoff_lrc
772 CALL pair_radius_setup(default_present, default_present, c_radius, c_radius, pair_radius)
774 subcells=subcells, nlname=
"sab_lrc")
775 CALL set_ks_env(ks_env=ks_env, sab_lrc=sab_lrc)
777 "/SAB_LRC",
"sab_lrc",
"SE LONG-RANGE CORRECTION")
783 IF (dft_control%qs_control%dftb_control%do_ewald)
THEN
784 CALL get_qs_env(qs_env=qs_env, ewald_env=ewald_env)
787 CALL pair_radius_setup(orb_present, orb_present, c_radius, c_radius, pair_radius)
789 subcells=subcells, nlname=
"sab_tbe")
790 CALL set_ks_env(ks_env=ks_env, sab_tbe=sab_tbe)
794 IF (dft_control%qs_control%dftb_control%dispersion)
THEN
795 IF (dft_control%qs_control%dftb_control%dispersion_type ==
dispersion_uff)
THEN
797 CALL get_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom)
800 default_present = .true.
801 CALL pair_radius_setup(default_present, default_present, c_radius, c_radius, pair_radius)
803 subcells=subcells, nlname=
"sab_vdw")
804 CALL set_ks_env(ks_env=ks_env, sab_vdw=sab_vdw)
809 IF (xtb .AND. (.NOT. dft_control%qs_control%xtb_control%do_tblite))
THEN
811 IF (dft_control%qs_control%xtb_control%do_ewald)
THEN
812 CALL get_qs_env(qs_env=qs_env, ewald_env=ewald_env)
815 CALL pair_radius_setup(orb_present, orb_present, c_radius, c_radius, pair_radius)
817 subcells=subcells, nlname=
"sab_tbe")
818 CALL set_ks_env(ks_env=ks_env, sab_tbe=sab_tbe)
821 pair_radius(1:nkind, 1:nkind) = dft_control%qs_control%xtb_control%rcpair(1:nkind, 1:nkind)
822 default_present = .true.
824 subcells=subcells, nlname=
"sab_xtb_pp")
825 CALL set_ks_env(ks_env=ks_env, sab_xtb_pp=sab_xtb_pp)
828 CALL get_qs_kind(qs_kind_set(ikind), xtb_parameter=xtb_atom)
831 default_present = .true.
832 CALL pair_radius_setup(default_present, default_present, c_radius, c_radius, pair_radius)
834 subcells=subcells, nlname=
"sab_xtbe")
835 CALL set_ks_env(ks_env=ks_env, sab_xtbe=sab_xtbe)
837 ALLOCATE (xb1_atom(nkind), xb2_atom(nkind))
838 c_radius = 0.5_dp*dft_control%qs_control%xtb_control%xb_radius
841 IF (zat == 17 .OR. zat == 35 .OR. zat == 53 .OR. zat == 85)
THEN
842 xb1_atom(ikind) = .true.
844 xb1_atom(ikind) = .false.
846 IF (zat == 7 .OR. zat == 8 .OR. zat == 15 .OR. zat == 16)
THEN
847 xb2_atom(ikind) = .true.
849 xb2_atom(ikind) = .false.
854 symmetric=.false., subcells=subcells, operator_type=
"PP", nlname=
"sab_xb")
857 "/SAB_XB",
"sab_xb",
"XB bonding")
860 IF (dft_control%qs_control%xtb_control%do_nonbonded &
861 .AND. (.NOT. dft_control%qs_control%xtb_control%do_tblite))
THEN
862 ngp =
SIZE(dft_control%qs_control%xtb_control%nonbonded%pot)
863 ALLOCATE (nonbond1_atom(nkind), nonbond2_atom(nkind))
864 nonbond1_atom = .false.
865 nonbond2_atom = .false.
868 rcut = sqrt(dft_control%qs_control%xtb_control%nonbonded%pot(ingp)%pot%rcutsq)
870 CALL get_atomic_kind(atomic_kind_set(ikind), element_symbol=element_symbol)
872 IF (trim(dft_control%qs_control%xtb_control%nonbonded%pot(ingp)%pot%at1) == trim(element_symbol))
THEN
873 nonbond1_atom(ikind) = .true.
875 CALL get_atomic_kind(atomic_kind_set(jkind), element_symbol=element_symbol2)
877 IF (trim(dft_control%qs_control%xtb_control%nonbonded%pot(ingp)%pot%at2) == trim(element_symbol2))
THEN
878 nonbond2_atom(jkind) = .true.
883 CALL pair_radius_setup(nonbond1_atom, nonbond2_atom, c_radius, c_radius, pair_radius)
885 symmetric=.false., subcells=subcells, operator_type=
"PP", nlname=
"sab_xtb_nonbond")
886 CALL set_ks_env(ks_env=ks_env, sab_xtb_nonbond=sab_xtb_nonbond)
887 CALL write_neighbor_lists(sab_xtb_nonbond, particle_set, cell, para_env, neighbor_list_section, &
888 "/SAB_XTB_NONBOND",
"sab_xtb_nonbond",
"XTB NONBONDED INTERACTIONS")
894 IF (.NOT. dft_control%qs_control%xtb_control%do_tblite)
THEN
895 CALL get_qs_env(qs_env=qs_env, dispersion_env=dispersion_env)
896 sab_vdw => dispersion_env%sab_vdw
897 sab_cn => dispersion_env%sab_cn
900 c_radius(:) = dispersion_env%rc_d4
902 c_radius(:) = dispersion_env%rc_disp
904 default_present = .true.
905 CALL pair_radius_setup(default_present, default_present, c_radius, c_radius, pair_radius)
907 subcells=subcells, operator_type=
"PP", nlname=
"sab_vdw")
908 dispersion_env%sab_vdw => sab_vdw
914 c_radius(ikind) = 4._dp*
ptable(zat)%covalent_radius*
bohr
916 CALL pair_radius_setup(default_present, default_present, c_radius, c_radius, pair_radius)
918 subcells=subcells, operator_type=
"PP", nlname=
"sab_cn")
919 dispersion_env%sab_cn => sab_cn
925 CALL get_qs_env(qs_env=qs_env, gcp_env=gcp_env)
926 IF (
ASSOCIATED(gcp_env))
THEN
927 IF (gcp_env%do_gcp)
THEN
928 sab_gcp => gcp_env%sab_gcp
930 c_radius(ikind) = gcp_env%gcp_kind(ikind)%rcsto
932 CALL pair_radius_setup(orb_present, orb_present, c_radius, c_radius, pair_radius)
934 subcells=subcells, operator_type=
"PP", nlname=
"sab_gcp")
935 gcp_env%sab_gcp => sab_gcp
937 NULLIFY (gcp_env%sab_gcp)
941 IF (lrigpw .OR. lri_optbas)
THEN
943 CALL pair_radius_setup(orb_present, orb_present, orb_radius, orb_radius, pair_radius)
944 soo_list => qs_env%lri_env%soo_list
946 mic=mic, molecular=molecule_only, subcells=subcells, nlname=
"soo_list")
947 qs_env%lri_env%soo_list => soo_list
949 "/SOO_LIST",
"soo_list",
"ORBITAL ORBITAL (RI)")
951 ALLOCATE (ri_present(nkind), ri_radius(nkind))
955 CALL get_qs_kind(qs_kind_set(ikind), basis_set=ri_basis_set, basis_type=
"RI_HXC")
956 IF (
ASSOCIATED(ri_basis_set))
THEN
957 ri_present(ikind) = .true.
960 ri_present(ikind) = .false.
964 CALL pair_radius_setup(orb_present, orb_present, orb_radius, orb_radius, pair_radius)
965 soo_list => qs_env%lri_env%soo_list
967 mic=mic, molecular=molecule_only, subcells=subcells, nlname=
"soo_list")
968 qs_env%lri_env%soo_list => soo_list
970 CALL pair_radius_setup(ri_present, ri_present, ri_radius, ri_radius, pair_radius)
971 saa_list => qs_env%lri_env%saa_list
973 mic=mic, molecular=molecule_only, subcells=subcells, nlname=
"saa_list")
974 qs_env%lri_env%saa_list => saa_list
976 CALL pair_radius_setup(ri_present, orb_present, ri_radius, orb_radius, pair_radius)
977 soa_list => qs_env%lri_env%soa_list
979 mic=mic, symmetric=.false., molecular=molecule_only, &
980 subcells=subcells, operator_type=
"ABC", nlname=
"saa_list")
981 qs_env%lri_env%soa_list => soa_list
987 CALL get_atomic_kind(atomic_kind_set(ikind), rcov=almo_rcov, rvdw=almo_rvdw)
989 c_radius(ikind) = max(almo_rcov, almo_rvdw)*
bohr* &
992 default_present = .true.
993 CALL pair_radius_setup(default_present, default_present, c_radius, c_radius, pair_radius)
995 subcells=subcells, operator_type=
"PP", nlname=
"sab_almo")
996 CALL set_ks_env(ks_env=ks_env, sab_almo=sab_almo)
1000 print_key_path =
"PRINT%DISTRIBUTION"
1005 basis_section=force_env_section, &
1006 print_key_path=print_key_path, &
1008 CALL write_neighbor_distribution(sab_orb, qs_kind_set, iw, para_env)
1011 basis_section=force_env_section, &
1012 print_key_path=print_key_path)
1019 DEALLOCATE (orb_present, default_present, core_present)
1020 DEALLOCATE (orb_radius, aux_fit_radius, c_radius, core_radius)
1021 DEALLOCATE (calpha, zeff)
1022 DEALLOCATE (pair_radius)
1023 IF (gth_potential_present .OR. sgp_potential_present)
THEN
1024 DEALLOCATE (ppl_present, ppl_radius)
1025 DEALLOCATE (ppnl_present, ppnl_radius)
1027 IF (paw_atom_present)
THEN
1028 DEALLOCATE (oce_present, oce_radius)
1030 IF (all_potential_present .OR. sgp_potential_present)
THEN
1031 DEALLOCATE (all_present, all_pot_rad)
1033 IF (cneo_potential_present)
THEN
1034 DEALLOCATE (cneo_present, nuc_orb_radius)
1037 CALL timestop(handle)
1065 mic, symmetric, molecular, subset_of_mol, current_subset, &
1066 operator_type, nlname, atomb_to_keep)
1073 REAL(
dp),
DIMENSION(:, :),
INTENT(IN) :: pair_radius
1074 REAL(
dp),
INTENT(IN) :: subcells
1075 LOGICAL,
INTENT(IN),
OPTIONAL :: mic, symmetric, molecular
1076 INTEGER,
DIMENSION(:),
OPTIONAL,
POINTER :: subset_of_mol
1077 INTEGER,
OPTIONAL :: current_subset
1078 CHARACTER(LEN=*),
INTENT(IN),
OPTIONAL :: operator_type
1079 CHARACTER(LEN=*),
INTENT(IN) :: nlname
1080 INTEGER,
DIMENSION(:),
INTENT(IN),
OPTIONAL :: atomb_to_keep
1082 CHARACTER(len=*),
PARAMETER :: routinen =
'build_neighbor_lists'
1085 INTEGER,
DIMENSION(:),
POINTER ::
list
1086 END TYPE local_lists
1088 INTEGER :: atom_a, atom_b, handle, i, iab, iatom, iatom_local, &
1089 iatom_subcell, icell, ikind, j, jatom, jatom_local, jcell, jkind, k, &
1090 kcell, maxat, mol_a, mol_b, nkind, otype, natom, inode, nnode, nentry
1091 INTEGER,
DIMENSION(3) :: cell_b, ncell, nsubcell, periodic
1092 INTEGER,
DIMENSION(:),
POINTER :: index_list
1093 LOGICAL :: include_ab, my_mic, &
1094 my_molecular, my_symmetric, my_sort_atomb
1095 LOGICAL,
ALLOCATABLE,
DIMENSION(:) :: pres_a, pres_b
1096 REAL(
dp) :: rab2, rab2_max, rab_max, rabm, deth, subcell_scale
1097 REAL(
dp),
DIMENSION(3) :: r, rab, ra, rb, sab_max, sb, &
1098 sb_pbc, sb_min, sb_max, rab_pbc, pd, sab_max_guard
1099 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: nlista, nlistb
1100 TYPE(local_lists),
DIMENSION(:),
POINTER :: lista, listb
1102 ALLOCATABLE,
DIMENSION(:) :: kind_a
1106 REAL(kind=
dp),
DIMENSION(:, :),
ALLOCATABLE :: r_pbc
1108 DIMENSION(:),
POINTER :: nl_iterator
1110 CALL timeset(routinen//
"_"//trim(nlname), handle)
1114 IF (
PRESENT(mic)) my_mic = mic
1115 my_symmetric = .true.
1116 IF (
PRESENT(symmetric)) my_symmetric = symmetric
1117 my_molecular = .false.
1119 IF (
PRESENT(molecular)) my_molecular = molecular
1121 IF (
PRESENT(operator_type))
THEN
1122 SELECT CASE (operator_type)
1127 cpassert(.NOT. my_molecular)
1128 my_symmetric = .false.
1131 my_symmetric = .false.
1141 my_sort_atomb = .false.
1142 IF (
PRESENT(atomb_to_keep))
THEN
1143 my_sort_atomb = .true.
1150 ALLOCATE (ab_list(nkind*nkind))
1151 DO iab = 1,
SIZE(ab_list)
1152 NULLIFY (ab_list(iab)%neighbor_list_set)
1153 ab_list(iab)%nl_size = -1
1154 ab_list(iab)%nl_start = -1
1155 ab_list(iab)%nl_end = -1
1156 NULLIFY (ab_list(iab)%nlist_task)
1160 ALLOCATE (pres_a(nkind), pres_b(nkind))
1162 pres_a(ikind) = any(pair_radius(ikind, :) > 0._dp)
1163 pres_b(ikind) = any(pair_radius(:, ikind) > 0._dp)
1167 natom =
SIZE(particle_set)
1168 ALLOCATE (r_pbc(3, natom))
1170 r_pbc(1:3, i) =
pbc(particle_set(i)%r(1:3), cell)
1176 maxat = max(maxat,
SIZE(
atom(ikind)%list))
1178 ALLOCATE (index_list(maxat))
1182 ALLOCATE (lista(nkind), listb(nkind), nlista(nkind), nlistb(nkind))
1186 NULLIFY (lista(ikind)%list, listb(ikind)%list)
1189 IF (
ASSOCIATED(
atom(ikind)%list_local_a_index))
THEN
1190 lista(ikind)%list =>
atom(ikind)%list_local_a_index
1191 nlista(ikind) =
SIZE(lista(ikind)%list)
1193 IF (
ASSOCIATED(
atom(ikind)%list_local_b_index))
THEN
1194 listb(ikind)%list =>
atom(ikind)%list_local_b_index
1195 nlistb(ikind) =
SIZE(listb(ikind)%list)
1198 IF (
ASSOCIATED(
atom(ikind)%list_local_a_index))
THEN
1199 lista(ikind)%list =>
atom(ikind)%list_local_a_index
1200 nlista(ikind) =
SIZE(lista(ikind)%list)
1202 nlistb(ikind) =
SIZE(
atom(ikind)%list)
1203 listb(ikind)%list => index_list
1205 CALL combine_lists(lista(ikind)%list, nlista(ikind), ikind,
atom)
1206 nlistb(ikind) =
SIZE(
atom(ikind)%list)
1207 listb(ikind)%list => index_list
1209 nlista(ikind) =
SIZE(
atom(ikind)%list_1d)
1210 lista(ikind)%list =>
atom(ikind)%list_1d
1211 nlistb(ikind) =
SIZE(
atom(ikind)%list)
1212 listb(ikind)%list => index_list
1221 maxat = max(maxat, nlista(ikind), nlistb(ikind))
1223 ALLOCATE (kind_a(2*maxat))
1226 CALL get_cell(cell=cell, periodic=periodic, deth=deth)
1230 IF (.NOT. pres_a(ikind)) cycle
1233 IF (.NOT. pres_b(jkind)) cycle
1235 iab = ikind + nkind*(jkind - 1)
1238 IF (pair_radius(ikind, jkind) <= 0._dp) cycle
1239 rab_max = pair_radius(ikind, jkind)
1240 IF (otype == 3)
THEN
1244 rabm = maxval(pair_radius(:, jkind))
1248 rab2_max = rabm*rabm
1255 sab_max_guard = 15.0_dp/pd
1258 subcell_scale = ((125.0_dp**3)/deth)**(1.0_dp/6.0_dp)
1262 nsubcell(:) = int(max(1.0_dp, min(0.5_dp*subcells*subcell_scale/sab_max(:), &
1263 0.5_dp*subcells*subcell_scale/sab_max_guard(:))))
1266 ncell(:) = (int(sab_max(:)) + 1)*periodic(:)
1269 symmetric=my_symmetric)
1270 neighbor_list_set => ab_list(iab)%neighbor_list_set
1272 DO iatom_local = 1, nlista(ikind)
1273 iatom = lista(ikind)%list(iatom_local)
1274 atom_a =
atom(ikind)%list(iatom)
1277 neighbor_list=kind_a(iatom_local)%neighbor_list)
1281 DO iatom_local = 1, nlista(ikind)
1282 iatom = lista(ikind)%list(iatom_local)
1283 atom_a =
atom(ikind)%list(iatom)
1284 r = r_pbc(:, atom_a)
1286 subcell(i, j, k)%natom = subcell(i, j, k)%natom + 1
1288 DO k = 1, nsubcell(3)
1289 DO j = 1, nsubcell(2)
1290 DO i = 1, nsubcell(1)
1291 maxat = subcell(i, j, k)%natom + subcell(i, j, k)%natom/10
1292 ALLOCATE (subcell(i, j, k)%atom_list(maxat))
1293 subcell(i, j, k)%natom = 0
1297 DO iatom_local = 1, nlista(ikind)
1298 iatom = lista(ikind)%list(iatom_local)
1299 atom_a =
atom(ikind)%list(iatom)
1300 r = r_pbc(:, atom_a)
1302 subcell(i, j, k)%natom = subcell(i, j, k)%natom + 1
1303 subcell(i, j, k)%atom_list(subcell(i, j, k)%natom) = iatom_local
1306 DO jatom_local = 1, nlistb(jkind)
1307 jatom = listb(jkind)%list(jatom_local)
1308 atom_b =
atom(jkind)%list(jatom)
1309 IF (my_sort_atomb .AND. .NOT. my_symmetric)
THEN
1310 IF (.NOT. any(atomb_to_keep == atom_b)) cycle
1312 IF (my_molecular)
THEN
1313 mol_b =
atom(jkind)%list_b_mol(jatom_local)
1314 IF (
PRESENT(subset_of_mol))
THEN
1315 IF (subset_of_mol(mol_b) .NE. current_subset) cycle
1318 r = r_pbc(:, atom_b)
1321 loop2_kcell:
DO kcell = -ncell(3), ncell(3)
1322 sb(3) = sb_pbc(3) + real(kcell,
dp)
1323 sb_min(3) = sb(3) - sab_max(3)
1324 sb_max(3) = sb(3) + sab_max(3)
1325 IF (periodic(3) /= 0)
THEN
1326 IF (sb_min(3) >= 0.5_dp)
EXIT loop2_kcell
1327 IF (sb_max(3) < -0.5_dp) cycle loop2_kcell
1331 loop2_jcell:
DO jcell = -ncell(2), ncell(2)
1332 sb(2) = sb_pbc(2) + real(jcell,
dp)
1333 sb_min(2) = sb(2) - sab_max(2)
1334 sb_max(2) = sb(2) + sab_max(2)
1335 IF (periodic(2) /= 0)
THEN
1336 IF (sb_min(2) >= 0.5_dp)
EXIT loop2_jcell
1337 IF (sb_max(2) < -0.5_dp) cycle loop2_jcell
1341 loop2_icell:
DO icell = -ncell(1), ncell(1)
1342 sb(1) = sb_pbc(1) + real(icell,
dp)
1343 sb_min(1) = sb(1) - sab_max(1)
1344 sb_max(1) = sb(1) + sab_max(1)
1345 IF (periodic(1) /= 0)
THEN
1346 IF (sb_min(1) >= 0.5_dp)
EXIT loop2_icell
1347 IF (sb_max(1) < -0.5_dp) cycle loop2_icell
1353 loop_k:
DO k = 1, nsubcell(3)
1354 loop_j:
DO j = 1, nsubcell(2)
1355 loop_i:
DO i = 1, nsubcell(1)
1359 IF (periodic(3) /= 0)
THEN
1360 IF (sb_max(3) < subcell(i, j, k)%s_min(3))
EXIT loop_k
1361 IF (sb_min(3) >= subcell(i, j, k)%s_max(3)) cycle loop_k
1364 IF (periodic(2) /= 0)
THEN
1365 IF (sb_max(2) < subcell(i, j, k)%s_min(2))
EXIT loop_j
1366 IF (sb_min(2) >= subcell(i, j, k)%s_max(2)) cycle loop_j
1369 IF (periodic(1) /= 0)
THEN
1370 IF (sb_max(1) < subcell(i, j, k)%s_min(1))
EXIT loop_i
1371 IF (sb_min(1) >= subcell(i, j, k)%s_max(1)) cycle loop_i
1374 IF (subcell(i, j, k)%natom == 0) cycle
1376 DO iatom_subcell = 1, subcell(i, j, k)%natom
1377 iatom_local = subcell(i, j, k)%atom_list(iatom_subcell)
1378 iatom = lista(ikind)%list(iatom_local)
1379 atom_a =
atom(ikind)%list(iatom)
1380 IF (my_molecular)
THEN
1381 mol_a =
atom(ikind)%list_a_mol(iatom_local)
1382 IF (mol_a /= mol_b) cycle
1384 IF (my_symmetric)
THEN
1385 IF (atom_a > atom_b)
THEN
1386 include_ab = (
modulo(atom_a + atom_b, 2) /= 0)
1388 include_ab = (
modulo(atom_a + atom_b, 2) == 0)
1390 IF (my_sort_atomb)
THEN
1391 IF ((.NOT. any(atomb_to_keep == atom_b)) .AND. &
1392 (.NOT. any(atomb_to_keep == atom_a)))
THEN
1393 include_ab = .false.
1399 IF (include_ab)
THEN
1400 ra(:) = r_pbc(:, atom_a)
1401 rab(:) = rb(:) - ra(:)
1402 rab2 = rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3)
1403 IF (rab2 < rab2_max)
THEN
1409 rab_pbc(:) =
pbc(rab(:), cell)
1410 IF (sum((rab_pbc - rab)**2) > epsilon(1.0_dp))
THEN
1411 include_ab = .false.
1414 IF (include_ab)
THEN
1416 neighbor_list=kind_a(iatom_local)%neighbor_list, &
1445 DEALLOCATE (lista(ikind)%list)
1450 DEALLOCATE (kind_a, pres_a, pres_b, lista, listb, nlista, nlistb)
1451 DEALLOCATE (index_list)
1458 IF (inode == 1) nentry = nentry + nnode
1462 ALLOCATE (ab_list(1)%nlist_task(nentry))
1463 ab_list(1)%nl_size = nentry
1464 DO iab = 2,
SIZE(ab_list)
1465 ab_list(iab)%nl_size = nentry
1466 ab_list(iab)%nlist_task => ab_list(1)%nlist_task
1475 iab = (ikind - 1)*nkind + jkind
1476 IF (ab_list(iab)%nl_start < 0) ab_list(iab)%nl_start = nentry
1477 IF (ab_list(iab)%nl_end < 0)
THEN
1478 ab_list(iab)%nl_end = nentry
1480 cpassert(ab_list(iab)%nl_end + 1 == nentry)
1481 ab_list(iab)%nl_end = nentry
1486 CALL timestop(handle)
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.