TABLE OF CONTENTS


ABINIT/m_double_grid [ Modules ]

[ Top ] [ Modules ]

NAME

  m_double_grid

FUNCTION

 This module defines the double grid object. This object contains the coarse mesh
 and the dense mesh used for the interpolation of the BSE Hamiltonian,
 and contains the mapping between the two meshes.

COPYRIGHT

 Copyright (C) 2008-2018 ABINIT group (YG, SP, MJV)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 MODULE m_double_grid
28 
29  use defs_basis
30  use m_errors
31  use m_abicore
32  use m_hide_blas
33  use m_bz_mesh
34  use m_kptrank
35 
36  use m_numeric_tools,  only : wrap2_zero_one, interpol3d_indices
37  use m_symtk,          only : matr3inv
38 
39  implicit none
40 
41  private

m_double_grid/compute_corresp [ Functions ]

[ Top ] [ m_double_grid ] [ Functions ]

NAME

 compute_corresp

FUNCTION

 Pre-process tables with mapping between divisions and coarse k-points

INPUTS

 double_grid

OUTPUT

 div2kdense(double_grid%nbz_coarse,double_grid%ndiv)
 (k_coarse,idiv) -> k_dense
 kdense2div(double_grid%nbz_dense)
 k_dense -> idiv

PARENTS

      m_hexc

CHILDREN

SOURCE

717 subroutine compute_corresp(double_grid, div2kdense, kdense2div)
718 
719 
720 !This section has been created automatically by the script Abilint (TD).
721 !Do not modify the following lines by hand.
722 #undef ABI_FUNC
723 #define ABI_FUNC 'compute_corresp'
724 !End of the abilint section
725 
726  implicit none
727 
728 !Argument ------------------------------------
729 !scalars
730  type(double_grid_t),intent(in) :: double_grid
731 !arrays
732  integer,intent(out) :: div2kdense(double_grid%nbz_coarse,double_grid%ndiv)
733  integer,intent(out) :: kdense2div(double_grid%nbz_dense)
734 
735 !Local variables -----------------------------
736 !scalars
737  integer :: iorder,ik_dense,ik_coarse
738 !arrays
739  integer :: curindices_dense(6)
740  integer,allocatable :: curindex(:)
741 
742 !*********************************************
743  ABI_MALLOC(curindex,(double_grid%nbz_coarse))
744  curindex = 1
745 
746  do ik_dense = 1,double_grid%nbz_dense
747 
748    ! From ik_ibz in the dense mesh -> indices_dense
749    iorder = double_grid%iktoint_dense(ik_dense)
750    !g01 = double_grid%g0_dense(:,iorder)
751 
752    ! From indices_dense -> indices_coarse
753    curindices_dense = double_grid%indices_dense(:,iorder)
754 
755    ik_coarse = double_grid%dense_to_coarse(ik_dense)
756    div2kdense(ik_coarse,curindex(ik_coarse)) = ik_dense
757    kdense2div(ik_dense) = curindex(ik_coarse)
758 
759    curindex(ik_coarse) = curindex(ik_coarse) + 1
760 
761  end do
762 
763  ABI_FREE(curindex)
764 
765 end subroutine compute_corresp

m_double_grid/compute_neighbours [ Functions ]

[ Top ] [ m_double_grid ] [ Functions ]

NAME

 compute_neighbours

FUNCTION

 Compute correspondance between points in the dense BZ and in the coarse BZ

INPUTS

   nbz_dense, nbz_closedcoarse, nbz_coarse, ndiv
   iktoint_dense(nbz_dense)
   indices_dense(6,nbz_dense)
   maxcomp_coarse(3)
   inttoik_coarse(nbz_closedcoarse)
   g0_coarse(3,nbz_closedcoarse)

OUTPUT

  dense_to_coarse(nbz_dense)
  coarse_to_dense(nbz_coarse,ndiv)

PARENTS

      m_double_grid

CHILDREN

SOURCE

630 subroutine compute_neighbours(nbz_dense, iktoint_dense, indices_dense, maxcomp_coarse, &
631 &  inttoik_coarse, g0_coarse, nbz_closedcoarse, nbz_coarse, ndiv, dense_to_coarse, coarse_to_dense)
632 
633 
634 !This section has been created automatically by the script Abilint (TD).
635 !Do not modify the following lines by hand.
636 #undef ABI_FUNC
637 #define ABI_FUNC 'compute_neighbours'
638 !End of the abilint section
639 
640  implicit none
641 
642 !Argument ------------------------------------
643 !scalars
644  integer,intent(in) :: nbz_dense, nbz_closedcoarse, nbz_coarse, ndiv
645 !arrays
646  integer,intent(in) :: iktoint_dense(nbz_dense)
647  integer,intent(in) :: indices_dense(6,nbz_dense)
648  integer,intent(in) :: maxcomp_coarse(3)
649  integer,intent(in) :: inttoik_coarse(nbz_closedcoarse)
650  integer,intent(in) :: g0_coarse(3,nbz_closedcoarse)
651  integer,intent(out) :: dense_to_coarse(nbz_dense)
652  integer,intent(out) :: coarse_to_dense(nbz_coarse,ndiv)
653 
654 !Local variables -----------------------------
655 !scalars
656  integer :: ik_dense, iorder, ik_coarse
657 !arrays
658  integer :: curindex(nbz_coarse)
659  integer :: curindices_dense(6), curindices_coarse(3)
660  integer :: g0(3)
661 
662 !*********************************************
663 
664  DBG_ENTER("COLL")
665 
666  coarse_to_dense = 1
667  dense_to_coarse = 1
668 
669  curindex = 1
670  do ik_dense = 1, nbz_dense
671   ! From ik_ibz in the dense mesh -> indices_dense
672   iorder = iktoint_dense(ik_dense)
673 
674   ! From indices_dense -> indices_coarse
675   curindices_dense = indices_dense(:,iorder)
676   curindices_coarse = curindices_dense(1:3)
677   ! From indices_coarse -> ik_ibz in the coarse mesh
678   call get_kpt_from_indices_coarse(curindices_coarse,maxcomp_coarse,&
679 &   inttoik_coarse,g0_coarse,nbz_closedcoarse,ik_coarse,g0)
680 
681   dense_to_coarse(ik_dense) = ik_coarse
682   coarse_to_dense(ik_coarse, curindex(ik_coarse)) = ik_dense
683 
684   curindex(ik_coarse) = curindex(ik_coarse) + 1
685  end do
686 
687  DBG_EXIT("COLL")
688 
689 end subroutine compute_neighbours

m_double_grid/create_indices_coarse [ Functions ]

[ Top ] [ m_double_grid ] [ Functions ]

NAME

 create_indices_coarse

FUNCTION

  Create mapping between kpoints and integer indexing

INPUTS

  bz(3,nbz) = k-points in the Brillouin Zone
  nbz = number of k-points
  klatt(3,3) = reciprocal space vectors defining the reciprocal cell
  nshiftk = Number of shifts
  shiftk(3,nshiftk) = Shiftks of the Brillouin Zone
  maxcomp(3) = Maximum int along each direction
  nbz_closed = Number of k-points inside the closed Brillouin Zone (adding periodic images)

OUTPUT

  indices(3,nbz_closed) = indices for each k-point in the closed BZ
  g0(3,nbz_closed) = g vectors between k-point inside bz and k-point given by indices
  iktoint(nbz) = mapping between k-points in the bz and int indices
  inttoik(nbz_closed) = mapping between int indices and k-points in the bz

PARENTS

      m_double_grid

CHILDREN

SOURCE

304 subroutine create_indices_coarse(bz, nbz, klatt, nshiftk, shiftk, maxcomp, nbz_closed, indices, g0, iktoint, inttoik)
305 
306 
307 !This section has been created automatically by the script Abilint (TD).
308 !Do not modify the following lines by hand.
309 #undef ABI_FUNC
310 #define ABI_FUNC 'create_indices_coarse'
311 !End of the abilint section
312 
313  implicit none
314 
315 !Argument ------------------------------------
316 !scalars
317  integer,intent(in) :: nbz,nshiftk,nbz_closed
318 !arrays
319  integer,intent(in) :: maxcomp(3)
320  integer,intent(out) :: indices(3,nbz_closed)
321  integer,intent(out) :: g0(3,nbz_closed)
322  integer,intent(out) :: iktoint(nbz), inttoik(nbz_closed)
323  real(dp),intent(in) :: bz(3,nbz),klatt(3,3),shiftk(3,nshiftk)
324 
325 !Local variables -----------------------------
326 !scalars
327  integer :: ik,ii,i1,i2, i3
328  logical :: found
329 !arrays
330  integer :: curg0(3)
331  real(dp) :: curk1(3),ktoget(3)
332 
333 !*********************************************
334 
335  ABI_CHECK(nshiftk==1,"nshiftk != 1 not supported")
336 
337  do i1 = 0,maxcomp(1)
338    do i2 = 0,maxcomp(2)
339      do i3 = 0,maxcomp(3)
340        ii = (i1*(maxcomp(2)+1)+i2)*(maxcomp(3)+1)+i3+1
341        ktoget(:) = shiftk(:,1)+(/i1,i2,i3/)
342        curk1(:) = MATMUL(klatt(:,:),ktoget(:))
343        found = .FALSE.
344        do ik = 1,nbz
345          if(isamek(curk1(:),bz(:,ik),curg0)) then
346            indices(:,ii) = (/i1,i2,i3/)
347            g0(:,ii) = curg0
348            if (i1 /= maxcomp(1) .and. i2 /= maxcomp(2) .and. i3 /= maxcomp(3)) then
349               iktoint(ik) = ii
350            end if
351            inttoik(ii) = ik
352            found = .TRUE.
353            exit
354          end if
355        end do
356        if (.not. found) then
357          write(std_out,*) "curk1 = ",curk1
358          write(std_out,*) bz
359          MSG_ERROR("A k-point generated from kptrlatt cannot be found in the BZ")
360        end if
361      end do
362    end do
363  end do
364 
365 end subroutine create_indices_coarse

m_double_grid/create_indices_dense [ Functions ]

[ Top ] [ m_double_grid ] [ Functions ]

NAME

 create_indices_dense

FUNCTION

  Create mapping between kpoints and integer indexing

INPUTS

  klatt_coarse(3,3) = reciprocal space vectors defining the reciprocal cell of coarse BZ
  maxcomp(3) = Maximum int along each direction
  bz_dense(3,nbz_dense) = k-points in the dense BZ
  nbz_dense = number of k-points in the dense BZ
  nshiftk = Number of shifts
  shiftk(3,nshiftk) = Shiftks of the Brillouin Zone
  kmult(3) = multiplication factors
  nbz_coarse = number of k-points in the coarse BZ
  kptrlatt_coarse(3,3) = real space vectors defining the reciprocal cell of coarse BZ

OUTPUT

  indices(6,nbz_dense) = indices for each k-point in the closed BZ
  g0(3,nbz_dense) = g vectors between k-point inside bz and k-point given by indices
  iktoint(nbz_dense) = mapping between k-points in the bz and int indices
  inttoik(nbz_dense) = mapping between int indices and k-points in the bz

PARENTS

      m_double_grid

CHILDREN

SOURCE

461 subroutine create_indices_dense(klatt_coarse, maxcomp, &
462 & bz_dense, nbz_dense, nshiftk, shiftk, kmult, indices, g0, inttoik, iktoint)
463 
464 
465 !This section has been created automatically by the script Abilint (TD).
466 !Do not modify the following lines by hand.
467 #undef ABI_FUNC
468 #define ABI_FUNC 'create_indices_dense'
469 !End of the abilint section
470 
471  implicit none
472 
473 !Argument ------------------------------------
474 !scalars
475  integer,intent(in) :: nbz_dense, nshiftk
476 !arrays
477  integer,intent(in) :: kmult(3),maxcomp(3)
478  integer,intent(out) :: indices(6,nbz_dense),g0(3,nbz_dense)
479  integer,intent(out) :: inttoik(nbz_dense),iktoint(nbz_dense)
480  real(dp),intent(in) :: bz_dense(3,nbz_dense),klatt_coarse(3,3)
481  real(dp),intent(in) :: shiftk(3,nshiftk)
482 
483 !Local variables -----------------------------
484  integer :: ik,ii,ii_coarse
485  integer :: i1,i2,i3,j1,j2,j3
486  logical :: found
487 !arrays
488  integer :: curg0(3)
489  real(dp) :: curk1(3),ktoget(3)
490 
491 !*********************************************
492 
493  call wrtout(std_out, "Create Indices Dense", "COLL")
494 
495  ABI_CHECK(nshiftk==1,"nshiftk != 1 not supported")
496 
497  do i1 = 0,maxcomp(1)-1
498    do i2 = 0,maxcomp(2)-1
499      do i3 = 0,maxcomp(3)-1
500        ii_coarse = (i1*(maxcomp(2)+1)+i2)*(maxcomp(3)+1)+i3+1
501 
502        do j1 = 0,kmult(1)-1
503          do j2 = 0,kmult(2)-1
504            do j3 = 0,kmult(3)-1
505              ii = ((i1*kmult(1)+j1)*(maxcomp(2)*kmult(2)) +&
506 &                  (i2*kmult(2)+j2))*(maxcomp(3)*kmult(3))+&
507 &                  (i3*kmult(3)+j3)+1
508 
509              ktoget(1) = i1+((REAL(j1)+shiftk(1,1))/kmult(1))
510              ktoget(2) = i2+((REAL(j2)+shiftk(2,1))/kmult(2))
511              ktoget(3) = i3+((REAL(j3)+shiftk(3,1))/kmult(3))
512 
513              curk1(:) = MATMUL(klatt_coarse(:,:),ktoget(:))
514              found = .FALSE.
515              do ik = 1,nbz_dense
516                if(isamek(curk1(:),bz_dense(:,ik),curg0)) then
517                  indices(:,ii) = (/i1,i2,i3,j1,j2,j3/)
518                  g0(:,ii) = curg0
519                  inttoik(ii) = ik
520                  iktoint(ik) = ii
521                  found = .TRUE.
522                  exit
523                end if
524              end do
525              if(.not. found) then
526                write(std_out,*) "curk1 = ",curk1
527                write(std_out,*) bz_dense
528                MSG_ERROR("Problem when creating indices")
529              end if
530            end do
531          end do
532        end do
533 
534      end do
535    end do
536  end do
537 
538 end subroutine create_indices_dense

m_double_grid/double_grid_free [ Functions ]

[ Top ] [ m_double_grid ] [ Functions ]

NAME

 double_grid_free

FUNCTION

 Deallocate all dynamics entities present in a double_grid structure.

INPUTS

 grid<double_grid>=The datatype to be freed.

SIDE EFFECTS

 All allocated memory is released.

PARENTS

      bethe_salpeter

CHILDREN

SOURCE

790 subroutine double_grid_free(grid)
791 
792 
793 !This section has been created automatically by the script Abilint (TD).
794 !Do not modify the following lines by hand.
795 #undef ABI_FUNC
796 #define ABI_FUNC 'double_grid_free'
797 !End of the abilint section
798 
799  implicit none
800 
801 !Arguments ------------------------------------
802  type(double_grid_t),intent(inout) :: grid
803 
804 ! *********************************************************************
805 
806  DBG_ENTER("COLL")
807 
808  !@double_grid
809 
810 !integer
811  if (allocated(grid%inttoik_coarse))  then
812    ABI_DEALLOCATE(grid%inttoik_coarse)
813  end if
814  if (allocated(grid%inttoik_dense))  then
815    ABI_DEALLOCATE(grid%inttoik_dense)
816  end if
817  if (allocated(grid%iktoint_coarse))  then
818    ABI_DEALLOCATE(grid%iktoint_coarse)
819  end if
820  if (allocated(grid%iktoint_dense))  then
821    ABI_DEALLOCATE(grid%iktoint_dense)
822  end if
823  if (allocated(grid%indices_coarse))  then
824    ABI_DEALLOCATE(grid%indices_coarse)
825  end if
826  if (allocated(grid%indices_dense))  then
827    ABI_DEALLOCATE(grid%indices_dense)
828  end if
829  if (allocated(grid%g0_coarse))  then
830    ABI_DEALLOCATE(grid%g0_coarse)
831  end if
832  if (allocated(grid%g0_dense))  then
833    ABI_DEALLOCATE(grid%g0_dense)
834  end if
835  if (allocated(grid%dense_to_coarse)) then
836    ABI_DEALLOCATE(grid%dense_to_coarse)
837  end if
838  if (allocated(grid%coarse_to_dense)) then
839    ABI_DEALLOCATE(grid%coarse_to_dense)
840  end if
841 
842 !real
843  if (allocated(grid%shiftk_dense)) then
844    ABI_DEALLOCATE(grid%shiftk_dense)
845  end if
846 
847  if (allocated(grid%shiftk_coarse)) then
848    ABI_DEALLOCATE(grid%shiftk_coarse)
849  end if
850 
851  DBG_EXIT("COLL")
852 
853 end subroutine double_grid_free

m_double_grid/double_grid_init [ Functions ]

[ Top ] [ m_double_grid ] [ Functions ]

NAME

 double_grid_init

FUNCTION

 Initialize the double_grid datatype "grid" from coarse and dense mesh

INPUTS

  Kmesh_coarse = descriptor of the coarse BZ sampling
  Kmesh_dense = descriptor of the dense BZ sampling
  kptrlatt_coarse(3,3) = vectors in R space that defines the reciprocal cell
  kmult(3) = multiplication factors from coarse to dense

OUTPUT

  grid = double_grid to be created

PARENTS

      setup_bse_interp

CHILDREN

SOURCE

177 subroutine double_grid_init(Kmesh_coarse,Kmesh_dense,kptrlatt_coarse,kmult,grid)
178 
179 
180 !This section has been created automatically by the script Abilint (TD).
181 !Do not modify the following lines by hand.
182 #undef ABI_FUNC
183 #define ABI_FUNC 'double_grid_init'
184 !End of the abilint section
185 
186  implicit none
187 
188 !Argument ------------------------------------
189 !scalars
190  type(kmesh_t),intent(in) :: Kmesh_coarse,Kmesh_dense
191  type(double_grid_t),intent(out) :: grid
192 !arrays
193  integer,intent(in) :: kptrlatt_coarse(3,3),kmult(3)
194 
195 !Local variables -----------------------------
196 !scalars
197  integer :: ii, info
198 !arrays
199  integer :: ipiv(3)
200  real(dp) :: rlatt_coarse(3,3),klatt_coarse(3,3),curmat(3,3)
201 
202 !*********************************************
203 
204  ABI_CHECK(Kmesh_coarse%nshift == 1, "Coarse mesh works only with nshiftk=1")
205  ABI_CHECK(Kmesh_dense%nshift == 1, "Dense mesh : Works only with nshiftk=1")
206 
207  grid%nshiftk_coarse = Kmesh_coarse%nshift
208  grid%nshiftk_dense = Kmesh_dense%nshift
209 
210  grid%nbz_coarse = Kmesh_coarse%nbz
211 
212  ABI_ALLOCATE(grid%shiftk_coarse,(3,grid%nshiftk_coarse))
213  ABI_ALLOCATE(grid%shiftk_dense,(3,grid%nshiftk_dense))
214 
215  grid%shiftk_coarse(:,:) = Kmesh_coarse%shift(:,:)
216  grid%shiftk_dense(:,:) = Kmesh_dense%shift(:,:)
217 
218  grid%kptrlatt_coarse(:,:) = kptrlatt_coarse(:,:)
219  rlatt_coarse(:,:) = kptrlatt_coarse(:,:)
220  call matr3inv(rlatt_coarse,klatt_coarse)
221  grid%klatt_coarse(:,:) = klatt_coarse(:,:)
222 
223  grid%nbz_dense = Kmesh_dense%nbz
224  grid%nbz_coarse = Kmesh_coarse%nbz
225 
226  grid%kmult(:) = kmult(:)
227  grid%ndiv = kmult(1)*kmult(2)*kmult(3)
228 
229  ABI_ALLOCATE(grid%indices_dense,(6,Kmesh_dense%nbz))
230  ABI_ALLOCATE(grid%g0_dense,(3,Kmesh_dense%nbz))
231  ABI_ALLOCATE(grid%iktoint_dense,(Kmesh_dense%nbz))
232  ABI_ALLOCATE(grid%inttoik_dense,(Kmesh_dense%nbz))
233 
234  grid%maxcomp_coarse(:) = -1
235 
236  curmat(:,:) = grid%kptrlatt_coarse(:,:)
237 
238  ! Gaussian elimination
239  call dgetrf(3,3,curmat,3,ipiv,info)
240 
241  grid%nbz_closedcoarse = 1
242 
243  do ii = 1,3
244    grid%maxcomp_coarse(ii) = ABS(NINT(curmat(ipiv(ii),ipiv(ii))))
245    grid%nbz_closedcoarse = grid%nbz_closedcoarse*(grid%maxcomp_coarse(ii)+1)
246  end do
247 
248  ABI_ALLOCATE(grid%indices_coarse,(3,grid%nbz_closedcoarse))
249  ABI_ALLOCATE(grid%g0_coarse,(3,grid%nbz_closedcoarse))
250  ABI_ALLOCATE(grid%iktoint_coarse,(Kmesh_coarse%nbz))
251  ABI_ALLOCATE(grid%inttoik_coarse,(grid%nbz_closedcoarse))
252 
253  ! We should pass 'grid' at this stage !
254 
255  call create_indices_coarse(Kmesh_coarse%bz, Kmesh_coarse%nbz, grid%klatt_coarse, &
256 &    grid%nshiftk_coarse, grid%shiftk_coarse, grid%maxcomp_coarse, grid%nbz_closedcoarse, grid%indices_coarse, &
257 &    grid%g0_coarse,grid%iktoint_coarse,grid%inttoik_coarse)
258 
259  call create_indices_dense(grid%klatt_coarse, grid%maxcomp_coarse, Kmesh_dense%bz, Kmesh_dense%nbz, &
260 &    grid%nshiftk_dense, grid%shiftk_dense, grid%kmult, grid%indices_dense, grid%g0_dense, grid%iktoint_dense, &
261 &    grid%inttoik_dense)
262 
263  ABI_ALLOCATE(grid%dense_to_coarse,(Kmesh_dense%nbz))
264  ABI_ALLOCATE(grid%coarse_to_dense,(Kmesh_coarse%nbz,grid%ndiv))
265 
266  call compute_neighbours(grid%nbz_dense, grid%iktoint_dense, grid%indices_dense, &
267 & grid%maxcomp_coarse, grid%inttoik_coarse, grid%g0_coarse, grid%nbz_closedcoarse, grid%nbz_coarse,&
268 & grid%ndiv, grid%dense_to_coarse, grid%coarse_to_dense)
269 
270 end subroutine double_grid_init

m_double_grid/double_grid_t [ Types ]

[ Top ] [ m_double_grid ] [ Types ]

NAME

 double_grid_t

FUNCTION

 The double grid contains a coarse mesh and a dense mesh
 It also contains the mapping between the two meshes

SOURCE

 54  type,public :: double_grid_t
 55 
 56   integer :: kmult(3)
 57   ! Number of subdivisions in the coarse box in each direction
 58 
 59   integer :: ndiv
 60   ! Total number of small sub-boxes in the large box associated to the coarse k-mesh.
 61   ! i.e. product(kmult)
 62 
 63   integer :: maxcomp_coarse(3)
 64   ! Dimensions of the box containing the coarse points in integer coord
 65 
 66   integer :: nbz_coarse
 67   ! Number of k-points in the coarse BZ (open mesh)
 68 
 69   integer :: nbz_closedcoarse
 70   ! Number of k-points inside the coarse BZ (closed mesh)
 71   ! = PROD(maxcomp_coarse+1)
 72 
 73   integer :: nbz_dense
 74   ! Number of k-point inside the dense BZ (open mesh)
 75   ! = PROD(maxcomp_coarse.*kmult)
 76 
 77   integer, allocatable :: inttoik_coarse(:)
 78   ! inttoik_coarse(nbz_closedcoarse)
 79   ! Index of the kpoint in the coarse BZ.
 80 
 81   integer, allocatable :: iktoint_coarse(:)
 82   ! iktoint_coarse(nbz_coarse)
 83   ! Index of int by kpoint
 84 
 85   integer, allocatable :: inttoik_dense(:)
 86 
 87   integer, allocatable :: iktoint_dense(:)
 88 
 89   integer,allocatable :: indices_coarse(:,:)
 90   ! indices_coarse(3,nbz_closedcoarse)
 91   ! Indices (i1,i2,i3) for each point, ordinated by the integer coord
 92 
 93   integer,allocatable :: indices_dense(:,:)
 94   ! indices_dense(6,nbz_dense)
 95   ! Indices (i1,i2,i3);(j1,j2,j3) for each point, ordinated by integer coord
 96 
 97   integer :: kptrlatt_dense(3,3)
 98   ! kptrlatt of the dense mesh
 99 
100   real(dp) :: klatt_dense(3,3)
101 
102   integer :: nshiftk_dense
103   ! Number of shifts in the dense mesh.
104 
105   real(dp),allocatable :: shiftk_dense(:,:)
106   ! shift_dense(3,nshiftk_dense)
107   ! Shifts of the dense mesh.
108 
109   ! Dense lattice
110   integer :: kptrlatt_coarse(3,3)
111   ! kptrlatt of the coarse mesh
112 
113   real(dp) :: klatt_coarse(3,3)
114 
115   integer :: nshiftk_coarse
116   ! Number of shifts in the coarse mesh.
117 
118   real(dp),allocatable :: shiftk_coarse(:,:)
119   ! shift_coarse(3,nshiftk_coarse)
120   ! Shifts of the coarse mesh.
121 
122   ! Coarse lattice
123   integer, allocatable :: g0_coarse(:,:)
124   integer, allocatable :: g0_dense(:,:)
125   ! g0_dense/coarse(3,nkpt_closedcoarse/dense)
126   ! G0 vector between the kpt obtained with indices
127   ! and the kpt obtained insize bz
128 
129   integer, allocatable :: dense_to_coarse(:)
130   ! dense_to_coarse(nbz_dense)
131   ! Give the ibz_coarse corresponding to the dense mesh (the (0,0,0) point)
132 
133   integer, allocatable :: coarse_to_dense(:,:)
134   ! coarse_to_dense(nbz_coarse,ndiv)
135   ! Give all the ibz_dense corresponding to the (0,0,0) coarse point
136 
137  end type double_grid_t
138 
139  public :: double_grid_init             ! Initializes the double grid with coarse mesh and dense mesh read from file
140  public :: double_grid_free             ! Deallocate all memory
141  public :: get_kpt_from_indices_coarse  ! Returns the k-point index and g0 vector associated to the set of indices
142  public :: compute_corresp              ! Compute correspondance data between k-dense and k-coarse
143  !public :: get_kpt_from_indices_dense
144 
145  public :: kptfine_av                   ! Find the k-points of a fine grid that are around a k-point of a coarse mesh.
146  public :: k_neighbors                  ! Find 8 neighbors of given k-point on a coarse grid, and return

m_double_grid/get_kpt_from_indices_coarse [ Functions ]

[ Top ] [ m_double_grid ] [ Functions ]

NAME

 get_kpt_from_indices_coarse

FUNCTION

  Returns the k-point index and g0 vector associated to the set of indices

INPUTS

  indices(3) = index of the searched k-point
  maxcomp(3) = Maximum int along each direction
  inttoik(nkpt) = mapping between int indices and k-points in the bz
  allg0(3,nkpt) = g vectors between k-point inside bz and k-point given by indices
  nkpt = number of k-points

OUTPUT

  ikpt = index of k-point we search
  g0(3) = g-vector obtained

PARENTS

      m_bseinterp,m_double_grid

CHILDREN

SOURCE

395 subroutine get_kpt_from_indices_coarse(indices,maxcomp,inttoik,allg0,nkpt,ikpt,g0)
396 !subroutine get_kpt_from_indices_coarse(double_grid,indices,ikpt,g0)
397 
398 !This section has been created automatically by the script Abilint (TD).
399 !Do not modify the following lines by hand.
400 #undef ABI_FUNC
401 #define ABI_FUNC 'get_kpt_from_indices_coarse'
402 !End of the abilint section
403 
404  implicit none
405 
406 !Argument ------------------------------------
407 !scalars
408  integer,intent(in) :: nkpt
409  integer,intent(out) :: ikpt
410 !arrays
411  integer,intent(in) :: indices(3),maxcomp(3)
412  integer,intent(in) :: inttoik(nkpt),allg0(3,nkpt)
413  integer,intent(out) :: g0(3)
414 
415 !Local variables -----------------------------
416 !scalars
417  integer :: curicoord
418 
419 !*********************************************
420 
421  curicoord = (indices(1)*(maxcomp(2)+1)+indices(2))*(maxcomp(3)+1)+indices(3)+1
422  ikpt = inttoik(curicoord)
423  g0 = allg0(:,curicoord)
424 
425 end subroutine get_kpt_from_indices_coarse

m_double_grid/get_kpt_from_indices_dense [ Functions ]

[ Top ] [ m_double_grid ] [ Functions ]

NAME

 get_kpt_from_indices_coarse

FUNCTION

  Returns the k-point index and g0 vector associated to the set of indices

INPUTS

  indices(6) = index of the searched k-point
  maxcomp(3) = Maximum int along each direction
  kmult(3) = multiplication factors
  inttoik(nkpt) = mapping between int indices and k-points in the bz
  allg0(3,nkpt) = g vectors between k-point inside bz and k-point given by indices
  nkpt = number of k-points

OUTPUT

  ikpt = index of k-point we search
  g0(3) = g-vector obtained

PARENTS

SOURCE

566 subroutine get_kpt_from_indices_dense(indices,maxcomp,kmult,inttoik,allg0,nkpt,ikpt,g0)
567 !subroutine get_kpt_from_indices_dense(double_grid,indices,ikpt,g0)
568 
569 !This section has been created automatically by the script Abilint (TD).
570 !Do not modify the following lines by hand.
571 #undef ABI_FUNC
572 #define ABI_FUNC 'get_kpt_from_indices_dense'
573 !End of the abilint section
574 
575  implicit none
576 
577 !Argument ------------------------------------
578 !scalars
579  integer, intent(in) :: nkpt
580  integer, intent(out) :: ikpt
581 !arrays
582  integer, intent(in) :: indices(6),maxcomp(3),inttoik(nkpt)
583  integer, intent(in) :: allg0(3,nkpt),kmult(3)
584  integer, intent(out) :: g0(3)
585 
586 !Local variables -----------------------------
587 !scalars
588  integer :: curicoord
589 
590 !*********************************************
591 
592  curicoord = ((indices(1)*kmult(1)+indices(4))*(maxcomp(2)*kmult(2))+&
593 &             (indices(2)*kmult(2)+indices(5)))*(maxcomp(3)*kmult(3))+&
594 &             (indices(3)*kmult(3)+indices(6))+1
595 
596  ikpt = inttoik(curicoord)
597  g0 = allg0(:,curicoord)
598 
599 end subroutine get_kpt_from_indices_dense

m_double_grid/k_neighbors [ Functions ]

[ Top ] [ m_double_grid ] [ Functions ]

NAME

   k_neighbors

FUNCTION

   find 8 neighbors of given k-point on a coarse grid, and return
   them along with relative k-shift within coarse grid cell

INPUTS

   kpt        = k-point to be interpolated to, in full BZ
   kptrlatt   = lattice vectors for coarse k-grid
   invrankkpt = rank list to find k-points

OUTPUT

   rel_kpt = k-point coordinates renormalized to coarse grid cell
   kpt_phon_indices = indices of k-points on corners of cell

TODO

  This routine is not used anymore. Deprecate or Remove?

PARENTS

CHILDREN

      get_rank_1kpt,interpol3d_indices,wrap2_zero_one

SOURCE

1050 subroutine k_neighbors (kpt, kptrlatt,kptrank_t, rel_kpt, kpt_phon_indices)
1051 
1052 
1053 !This section has been created automatically by the script Abilint (TD).
1054 !Do not modify the following lines by hand.
1055 #undef ABI_FUNC
1056 #define ABI_FUNC 'k_neighbors'
1057 !End of the abilint section
1058 
1059  implicit none
1060 
1061 ! inputs
1062  real(dp), intent(in) :: kpt(3)
1063  integer, intent(in) :: kptrlatt(3,3)
1064  type(kptrank_type), intent(in) :: kptrank_t
1065 
1066 ! outputs
1067  real(dp), intent(out) :: rel_kpt(3)
1068  integer, intent(out) :: kpt_phon_indices(8)
1069 
1070 ! local vars
1071  integer :: symrankkpt
1072  integer :: ir1,ir2,ir3, pr1,pr2,pr3
1073  real(dp) :: redkpt(3), cornerkpt(3), res
1074 
1075 ! *************************************************************************
1076 
1077 !wrap fine kpt to [0,1]
1078  call wrap2_zero_one(kpt(1),redkpt(1),res)
1079  call wrap2_zero_one(kpt(2),redkpt(2),res)
1080  call wrap2_zero_one(kpt(3),redkpt(3),res)
1081 !find 8 indices of points neighboring ikpt_phon, for interpolation
1082  call interpol3d_indices (redkpt,kptrlatt(1,1),kptrlatt(2,2),kptrlatt(3,3), &
1083 & ir1,ir2,ir3, pr1,pr2,pr3)
1084 
1085 !transpose ir pr to ikpt_phon indices
1086 !order of kpt_phons:
1087 !ir1 ir2 ir3
1088  cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/)
1089  call get_rank_1kpt (cornerkpt,symrankkpt,kptrank_t)
1090  kpt_phon_indices(1) = kptrank_t%invrank(symrankkpt)
1091 !pr1 ir2 ir3
1092  cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/)
1093  call get_rank_1kpt (cornerkpt,symrankkpt,kptrank_t)
1094  kpt_phon_indices(2) = kptrank_t%invrank(symrankkpt)
1095 !ir1 pr2 ir3
1096  cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/)
1097  call get_rank_1kpt (cornerkpt,symrankkpt,kptrank_t)
1098  kpt_phon_indices(3) = kptrank_t%invrank(symrankkpt)
1099 !pr1 pr2 ir3
1100  cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/)
1101  call get_rank_1kpt (cornerkpt,symrankkpt,kptrank_t)
1102  kpt_phon_indices(4) = kptrank_t%invrank(symrankkpt)
1103 !ir1 ir2 pr3
1104  cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/)
1105  call get_rank_1kpt (cornerkpt,symrankkpt,kptrank_t)
1106  kpt_phon_indices(5) = kptrank_t%invrank(symrankkpt)
1107 !pr1 ir2 pr3
1108  cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/)
1109  call get_rank_1kpt (cornerkpt,symrankkpt,kptrank_t)
1110  kpt_phon_indices(6) = kptrank_t%invrank(symrankkpt)
1111 !ir1 pr2 pr3
1112  cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/)
1113  call get_rank_1kpt (cornerkpt,symrankkpt,kptrank_t)
1114  kpt_phon_indices(7) = kptrank_t%invrank(symrankkpt)
1115 !pr1 pr2 pr3
1116  cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/)
1117  call get_rank_1kpt (cornerkpt,symrankkpt,kptrank_t)
1118  kpt_phon_indices(8) = kptrank_t%invrank(symrankkpt)
1119 
1120 !retrieve the gkq matrix for all q, at the neighbor k vectors
1121  rel_kpt(1) = redkpt(1)*kptrlatt(1,1)-real(ir1-1)
1122  rel_kpt(2) = redkpt(2)*kptrlatt(2,2)-real(ir2-1)
1123  rel_kpt(3) = redkpt(3)*kptrlatt(3,3)-real(ir3-1)
1124 
1125 end subroutine k_neighbors

m_double_grid/kptfine_av [ Functions ]

[ Top ] [ m_double_grid ] [ Functions ]

NAME

 kptfine_av

FUNCTION

 Find the k-points of a fine grid that are around a k-point of a coarse mesh.

INPUTS

  center(3) = the point of the coarse mesh around which you want know which
              k-points of the fine mesh belong to.
  qptrlatt(3,3) = qptrlatt of the considered calculation (this is obtained
              from the input variable ngqpt and shiftq.
  kpt_fine(3,nkpt_fine) = this table contain all the k-points of the fine grid
              in the full BZ (no sym op. allowed) and is read from the header
              of the dense WF file.
  nkpt_fine = number of k-points of the fine grid read from the header of the
              dense WF file.

OUTPUT

  kpt_fine_sub(nkpt_sub) = k-points of the fine grid that are around center(3)
  nkpt_sub = number of k-points of the fine grid that are around center(3)
  wgt_sub(nkpt_sub) = weight of the k-points of the fine grid that are around center(3).

PARENTS

      eig2stern,eig2tot

CHILDREN

SOURCE

 888 subroutine kptfine_av(center,qptrlatt,kpt_fine,nkpt_fine,kpt_fine_sub,nkpt_sub,wgt_sub)
 889 
 890 
 891 !This section has been created automatically by the script Abilint (TD).
 892 !Do not modify the following lines by hand.
 893 #undef ABI_FUNC
 894 #define ABI_FUNC 'kptfine_av'
 895 !End of the abilint section
 896 
 897  implicit none
 898 
 899 !Arguments ------------------------------------
 900 !scalars
 901  integer,intent(in)   :: nkpt_fine
 902  integer,intent(out)  :: nkpt_sub
 903 !arrays
 904  integer,intent(in)   :: qptrlatt(3,3)
 905  real(dp),intent(in)  :: kpt_fine(3,nkpt_fine)
 906  real(dp),intent(in)  :: center(3)
 907  integer,pointer      :: kpt_fine_sub(:)
 908  real(dp),pointer     :: wgt_sub(:)
 909 
 910 !Local variables-------------------------------
 911 !scalars
 912  integer :: ikpt,aa,bb,cc
 913  integer :: ii,jj
 914 !arrays
 915  real(dp) :: center_ref(3)
 916  real(dp) :: kpt_fine_ref(3)
 917  real(dp) :: kpt_tmp(3),kpt_tmp2(3)
 918  integer,allocatable  :: kpt_fine_sub_tmp(:)
 919  real(dp),allocatable :: wgt_sub_tmp(:)
 920  logical :: found(3)
 921 
 922 ! *************************************************************************
 923 
 924  ABI_ALLOCATE(kpt_fine_sub_tmp,(nkpt_fine))
 925  ABI_ALLOCATE(wgt_sub_tmp,(nkpt_fine))
 926 
 927 !It is easier to work in real space using the qptrlatt matrices because in this
 928 !referential any k-points sampling will be cast into an orthorhombic shape.
 929 !In that space we can simply take all k-points of the fine grid that between
 930 !center_ref-0.5 and center_ref+0.5
 931 
 932  center_ref = MATMUL(qptrlatt,center)
 933 
 934 !When considering points center(3) that lying close or on a BZ edge we need to
 935 !take the k-points of the fine grid taking into account unklamp vectors. This
 936 !is done with the aa, bb and cc loops.
 937 
 938  ii = 1
 939  do ikpt=1,nkpt_fine
 940    kpt_tmp = kpt_fine(:,ikpt)
 941    do aa=-1,1
 942      kpt_tmp2(1) = kpt_tmp(1)+aa
 943      do bb=-1,1
 944        kpt_tmp2(2) = kpt_tmp(2)+bb
 945        do cc=-1,1
 946          kpt_tmp2(3) = kpt_tmp(3)+cc
 947          kpt_fine_ref = MATMUL(qptrlatt,kpt_tmp2)
 948          if((kpt_fine_ref(1)>=center_ref(1)-0.5-tol8).and.&
 949 &         (kpt_fine_ref(1)<=center_ref(1)+0.5+tol8)) then
 950            if((kpt_fine_ref(2)>=center_ref(2)-0.5-tol8).and.&
 951 &           (kpt_fine_ref(2)<=center_ref(2)+0.5+tol8)) then
 952              if((kpt_fine_ref(3)>=center_ref(3)-0.5-tol8).and.&
 953 &             (kpt_fine_ref(3)<=center_ref(3)+0.5+tol8)) then
 954                kpt_fine_sub_tmp(ii) = ikpt
 955                ii = ii +1
 956              end if
 957            end if
 958          end if
 959        end do
 960      end do
 961    end do
 962  end do
 963 
 964  nkpt_sub = ii-1
 965  ABI_ALLOCATE(kpt_fine_sub,(nkpt_sub))
 966  ABI_ALLOCATE(wgt_sub,(nkpt_sub))
 967 
 968  do jj=1,nkpt_sub
 969    kpt_fine_sub(jj) = kpt_fine_sub_tmp(jj)
 970  end do
 971 
 972 !We then compute a weight function. This weight function is simply a
 973 !rectangular weight function that take the value 1 for k-points of the fine
 974 !grid inside the cube, 0.5 for k-points that are lying on one face of the cube,
 975 !0.25 for k-points that are lying on an edge of the cube and 0.125 for k-points
 976 !that are lying on a peak of the cube.
 977 
 978  wgt_sub(:) = 1.0
 979 
 980  do ikpt=1,nkpt_sub
 981    found(:) = .True.
 982    kpt_tmp = kpt_fine(:,kpt_fine_sub(ikpt))
 983    do aa=-1,1
 984      kpt_tmp2(1) = kpt_tmp(1)+aa
 985      do bb=-1,1
 986        kpt_tmp2(2) = kpt_tmp(2)+bb
 987        do cc=-1,1
 988          kpt_tmp2(3) = kpt_tmp(3)+cc
 989          kpt_fine_ref = MATMUL(qptrlatt,kpt_tmp2)
 990          if((ABS(kpt_fine_ref(1)-center_ref(1)-0.5)< tol8) .or.&
 991 &         (ABS(kpt_fine_ref(1)-center_ref(1)+0.5) < tol8)) then
 992            if(found(1)) then
 993              wgt_sub(ikpt) = wgt_sub(ikpt)*0.5
 994              found(1) = .False.
 995            end if
 996          end if
 997          if((ABS(kpt_fine_ref(2)-center_ref(2)-0.5) < tol8) .or.&
 998 &         (ABS(kpt_fine_ref(2)-center_ref(2)+0.5) < tol8)) then
 999            if(found(2)) then
1000              wgt_sub(ikpt) = wgt_sub(ikpt)*0.5
1001              found(2) = .False.
1002            end if
1003          end if
1004          if((ABS(kpt_fine_ref(3)-center_ref(3)-0.5)< tol8) .or.&
1005 &         (ABS(kpt_fine_ref(3)-center_ref(3)+0.5) < tol8)) then
1006            if(found(3)) then
1007              wgt_sub(ikpt) = wgt_sub(ikpt)*0.5
1008              found(3) = .False.
1009            end if
1010          end if
1011        end do
1012      end do
1013    end do
1014  end do
1015 
1016  ABI_DEALLOCATE(kpt_fine_sub_tmp)
1017  ABI_DEALLOCATE(wgt_sub_tmp)
1018 
1019 end subroutine kptfine_av