xy_shells Subroutine

public pure subroutine xy_shells(lattice, shell_table, distance_set)

Build a /ordered/ set of all inter-atomic distances in a given xy_lattice and a table storing distances among all sites. Searching this matrix for the n-th set entry gives a mask of all pairs of atomic sites that are n-th order neigbors. E.g. shell_table == distance_set(1) would provide a mask for nearest neighbors, from which we can build the tight-binding hopping term of the corresponding lattice hamiltonian.

Arguments

Type IntentOptional Attributes Name
type(xy_lattice), intent(in) :: lattice
real(kind=8), intent(out), allocatable :: shell_table(:,:)
real(kind=8), intent(out), allocatable :: distance_set(:)

Calls

proc~~xy_shells~~CallsGraph proc~xy_shells xy_neighbors::xy_shells proc~xy_distance xy_coordinates::xy_distance proc~xy_shells->proc~xy_distance sort sort proc~xy_shells->sort proc~xy_norm xy_coordinates::xy_norm proc~xy_distance->proc~xy_norm

Called by

proc~~xy_shells~~CalledByGraph proc~xy_shells xy_neighbors::xy_shells proc~xy_nearest_neighbors xy_neighbors::xy_nearest_neighbors proc~xy_nearest_neighbors->proc~xy_shells proc~xy_next_nearest_neighbors xy_neighbors::xy_next_nearest_neighbors proc~xy_next_nearest_neighbors->proc~xy_shells

Source Code

   pure subroutine xy_shells(lattice,shell_table,distance_set)
      !! Build a /ordered/ set of all inter-atomic distances
      !! in a given xy_lattice and a table storing distances
      !! among all sites. Searching this matrix for the n-th
      !! set entry gives a mask of all pairs of atomic sites
      !! that are n-th order neigbors.
      !! E.g. shell_table == distance_set(1) would provide a
      !! mask for nearest neighbors, from which we can build
      !! the tight-binding hopping term of the corresponding
      !! lattice hamiltonian.
      type(xy_lattice),intent(in)      :: lattice
      real(8),allocatable,intent(out)  :: distance_set(:)
      real(8),allocatable,intent(out)  :: shell_table(:,:)
      real(8)                          :: distance
      type(xy_site)                    :: A,B
      integer                          :: L,i,j,k
      L = size(lattice%site)
      allocate(shell_table(L,L),distance_set(L**2))
      shell_table = 0d0 ! init the shells
      distance_set = 0d0 ! and the set
      k = 0 ! unique distance counter
      do i = 1,L-1
         do j = i+1,L
            A = lattice%site(i)
            B = lattice%site(j)
            distance = xy_distance(A,B)
            shell_table(i,j) = distance
            shell_table(j,i) = distance
            if(all(abs(distance_set-distance)>1d-12))then
               k = k + 1
               distance_set(k) = distance
            end if
         enddo
      enddo
      ! Shrink the set to the actual required size:
      distance_set = pack(distance_set,distance_set/=0d0)
      call sort(distance_set)
   end subroutine xy_shells