1module m_supercell 2 3! Calculates the supercell factors by examining all atoms. 4 5 use precision, only : dp 6 7 implicit none 8 9 public :: exact_sc_ag 10 private 11 12contains 13 14 !> Calculate the exact super-cell size based on the 15 !> atomic connection graph. The logic is similar to that in [[hsparse]], 16 !> but instead of orbitals we deal with atoms. 17 !> 18 !> Two atoms are connected if any of their orbitals overlap, including 19 !> the possibility of third-center KB projectors mediating an indirect overlap. 20 !> For a pictorial representation of this graph, 21 !> see [this page](|page|/implementation/1-auxiliary-supercell.html) 22 !> 23 !> Since this routine calculates explicitly the atom connections, 24 !> the resulting super-cell size will be exact. 25 26 subroutine exact_sc_ag(negl,ucell,na_u,isa,xa,nsc) 27 28 use class_iSpData2D 29 30 use atom_graph 31 32#ifdef MPI 33 use class_OrbitalDistribution 34 use parallel, only : Nodes 35 use parallelsubs, only : set_blocksizedefault 36 use mpi_siesta 37#endif 38 39 !> Whether the KB projectors are neglected or not 40 logical, intent(in) :: negl 41 !> Unit-cell vectors 42 real(dp), intent(in) :: ucell(3,3) 43 !> Total number of atoms in the unit cell 44 integer, intent(in) :: na_u 45 !> Atomic species array 46 integer, intent(in) :: isa(na_u) 47 !> Atomic coordinates 48 real(dp), intent(in) :: xa(3,na_u) 49 !> Supercell multipliers along the three unit-cell vectors 50 integer, intent(out) :: nsc(3) 51 52 ! Calculate the atomic graph 53 type(tAtomGraph) :: ag 54#ifdef MPI 55 type(OrbitalDistribution) :: dit 56#endif 57 integer, pointer :: sc(:,:) 58 integer :: ia, na_l, n_nzs 59 60#ifdef MPI 61 if ( na_u > Nodes ) then ! We can only distribute something if all have something 62 call set_blocksizedefault(Nodes,na_u,ia) 63 call newDistribution(ia,MPI_Comm_World,dit,name='AG-dist') 64 na_l = num_local_elements(dit,na_u) 65 call atom_graph_generate( negl, ucell, na_u, isa, xa, ag, dit , & 66 set_xijo = .false. ) 67 else 68 call atom_graph_generate( negl, ucell, na_u, isa, xa, ag , & 69 set_xijo = .false. ) 70 end if 71#else 72 na_l = na_u 73 call atom_graph_generate( negl, ucell, na_u, isa, xa, ag , & 74 set_xijo = .false. ) 75#endif 76 77 ! Extract the largest super-cell index from the sc variable. 78 n_nzs = nnzs(ag%sc_2D) 79 sc => val(ag%sc_2D) 80 81 !> Find the biggest supercell 82 !> We make it symmetric on both sides. 83 nsc(:) = 0 84 do ia = 1 , n_nzs 85 if ( abs(sc(1,ia)) > nsc(1) ) nsc(1) = abs(sc(1,ia)) 86 if ( abs(sc(2,ia)) > nsc(2) ) nsc(2) = abs(sc(2,ia)) 87 if ( abs(sc(3,ia)) > nsc(3) ) nsc(3) = abs(sc(3,ia)) 88 end do 89 90 ! DEBUG 91 !call atom_graph_print(ag,na_u,isa,xa) 92 93 ! Clean-up 94 call delete(ag) 95 96#ifdef MPI 97 ! Reduce the nsc 98 nullify(sc) ; allocate(sc(3,1)) 99 call MPI_AllReduce(nsc(1),sc(1,1),3,MPI_Integer,MPI_MAX, & 100 MPI_Comm_World, ia) 101 nsc(:) = sc(:,1) 102 deallocate(sc) 103 call delete(dit) 104#endif 105 106 ! nsc now contains the one-sided range of image-cell neighbors in each of the 107 ! unit-cell directions. 108 ! Hence, the correct number of super-cells is: 109 nsc(:) = nsc(:) * 2 + 1 110 111 end subroutine exact_sc_ag 112 113end module m_supercell 114