module params
   implicit none
   ! x and y dimensions of the lattice, Lx=even
   integer,parameter::Lx=8,Ly=4, LL=Lx*Ly
   integer,dimension(0:LL-1)::ann,bnn,cnn
end module

program graphene_neighbors
   use params
   implicit none
   integer::j
   call tabulate_neighbors()
   write(*,"(4a10)"),"j","ann(j)","bnn(j)","cnn(j)"; print*
   do j=0,LL-1
      write(*,"(4i10)"),j,ann(j),bnn(j),cnn(j)
   enddo
end program

subroutine tabulate_neighbors()
   use params
   implicit none
   integer::j
   ! this subroutine returns n.n.'s in a,b, and c directions
   ! see file honeycomb-lattice.png 
   do j=0,LL-2,2        ! even sites
      ann(j)=j-1
      bnn(j)=j+1
      cnn(j)=mod(j+Lx-1,LL)
   enddo
   ann(0)=LL-1
   do j=1,LL-1,2        ! odd sites
      ann(j)=j+1
      bnn(j)=j-1
      cnn(j)=mod(j-(Lx-1)+LL,LL)
   enddo
   ann(LL-1)=0
   
end subroutine
