TABLE OF CONTENTS


ABINIT/m_bseinterp [ Modules ]

[ Top ] [ Modules ]

NAME

 m_bseinterp

FUNCTION

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (M.Giantomassi, Y. Gillet)
  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

CHILDREN

SOURCE

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 #include "abi_common.h"
25 
26 MODULE m_bseinterp
27 
28  use defs_basis
29  use m_abicore
30  use m_bs_defs
31  use m_xmpi
32  use m_errors
33  use m_nctk
34  use m_haydock_io
35  use m_linalg_interfaces
36 #ifdef HAVE_NETCDF
37  use netcdf
38 #endif
39 
40  use m_fstrings,          only : indent, strcat, sjoin, itoa
41  use defs_datatypes,      only : pseudopotential_type
42  use m_hide_blas,         only : xdotc
43  use m_fft_mesh,          only : calc_ceigr
44  use m_crystal,           only : crystal_t
45  use m_bz_mesh,           only : kmesh_t
46  use m_double_grid,       only : double_grid_t, get_kpt_from_indices_coarse
47  use m_wfd,               only : wfd_t, wfd_sym_ur, wfd_get_ur, wfd_change_ngfft
48  use m_pawtab,            only : pawtab_type
49 
50  implicit none
51 
52  private

m_bseinterp/int_alloc_work [ Functions ]

[ Top ] [ m_bseinterp ] [ Functions ]

NAME

 int_alloc_work

FUNCTION

 Allocate temporary arrays

INPUTS

OUTPUT

PARENTS

      m_hexc

CHILDREN

SOURCE

242 subroutine int_alloc_work(interpolator, work_size)
243 
244 
245 !This section has been created automatically by the script Abilint (TD).
246 !Do not modify the following lines by hand.
247 #undef ABI_FUNC
248 #define ABI_FUNC 'int_alloc_work'
249 !End of the abilint section
250 
251  implicit none
252 
253 !Arguments ---------------------------
254 !scalars
255  integer,intent(in) :: work_size
256  type(interpolator_t),intent(inout) :: interpolator
257 
258 !*****************************************************************************
259 
260  ABI_MALLOC(interpolator%btemp,(work_size))
261  ABI_MALLOC(interpolator%ctemp,(work_size))
262 
263 end subroutine int_alloc_work

m_bseinterp/int_compute_overlaps [ Functions ]

[ Top ] [ m_bseinterp ] [ Functions ]

NAME

 int_compute_overlaps

FUNCTION

 Compute the overlaps prefactors

INPUTS

OUTPUT

PARENTS

      m_bseinterp

CHILDREN

SOURCE

333 subroutine int_compute_overlaps(interpolator, double_grid, Wfd_dense, Wfd_coarse, &
334 &   Kmesh_dense, Kmesh_coarse, BSp, Cryst, Psps, Pawtab)
335 
336 
337 !This section has been created automatically by the script Abilint (TD).
338 !Do not modify the following lines by hand.
339 #undef ABI_FUNC
340 #define ABI_FUNC 'int_compute_overlaps'
341 !End of the abilint section
342 
343  implicit none
344 
345 !Arguments ---------------------------
346 !scalars
347  type(interpolator_t),intent(inout) :: interpolator
348  type(double_grid_t),intent(in),target :: double_grid
349  type(wfd_t),intent(inout) :: Wfd_dense, Wfd_coarse
350  type(kmesh_t),intent(in) :: Kmesh_dense, Kmesh_coarse
351  type(excparam),intent(in) :: BSp
352  type(crystal_t),intent(in) :: Cryst
353  type(pseudopotential_type),intent(in) :: Psps
354 !arrays
355  type(pawtab_type),intent(in) :: Pawtab(Cryst%ntypat*Wfd_coarse%usepaw)
356 
357 !Local variables ---------------------
358 !scalars
359  integer :: nprocs, my_rank
360  integer :: ierr
361  integer :: nfft, nspinor, nsppol, nvert
362  integer :: ib_coarse, ib_dense
363  integer :: ik_coarse, ik_dense
364  integer :: spin
365  integer :: iorder
366  integer :: ivertex, ix, iy, iz
367  integer :: bstart, bstop
368  real(dp),parameter :: threshold = 0.1_dp
369  complex(gwpc) :: ovlp
370 !arrays
371  integer :: curindices_dense(6), curindices_coarse(3)
372  integer :: neighbour(3)
373  integer :: g0(3),g01(3),diffg0(3)
374  complex(gwpc),allocatable :: ur_coarse(:),ur_dense(:)
375  complex(gwpc),allocatable :: ceigr(:)
376 !arrays
377 
378 !*****************************************************************************
379 
380  nprocs = xmpi_comm_size(Wfd_coarse%comm)
381  my_rank = xmpi_comm_rank(Wfd_coarse%comm)
382 
383  ABI_UNUSED(Pawtab(1)%basis_size)
384 
385  ! Ensure Wfd and Wfd_coarse use the same FFT mesh.
386  call wfd_change_ngfft(Wfd_dense,Cryst,Psps,Wfd_coarse%ngfft)
387  nfft = Wfd_coarse%nfft
388  nspinor = Wfd_coarse%nspinor
389  nsppol = Bsp%nsppol
390  nvert = interpolator%nvert
391 
392  ! Allocate workspace for wavefunctions in real space.
393  ABI_MALLOC(ur_coarse,(nfft*nspinor))
394  ABI_MALLOC(ur_dense,(nfft*nspinor))
395  ABI_MALLOC(ceigr,(nfft*nspinor))
396 
397  interpolator%overlaps = czero
398 
399  ! TODO
400  ! 1) Choose whether we want to compute only dvv, dcc or all dbb
401  ! 2) Check the ordering of the loops
402  ! 3) Improve vertex -> neighbour (in double_grid ?)
403  do spin = 1,nsppol
404    do ik_dense = 1,double_grid%nbz_dense
405 
406      ! MPI parallelization
407      ! We assume that each node owns in memory the full set of wavefunctions
408      ! both coarse and dense k-mesh and both spins.
409      if (mod(ik_dense, nprocs) /= my_rank) cycle
410 
411      ! From ik_dense -> indices_dense
412      iorder = double_grid%iktoint_dense(ik_dense)
413      g01 = double_grid%g0_dense(:,iorder)
414      curindices_dense = double_grid%indices_dense(:,iorder)
415 
416      do ivertex = 1,nvert
417 
418        ! From vertex to neighbour
419        ! TODO improve this part + permit to choose other neighbour (e.g. nearest neighbour for RL)
420        if(nvert > 1) then
421          ix = (ivertex-1)/4
422          iy = (ivertex-ix*4-1)/2
423          iz = (ivertex-ix*4-iy*2-1)
424        else
425          ix = (BSp%rl_nb-1)/4
426          iy = (BSp%rl_nb-ix*4-1)/2
427          iz = (BSp%rl_nb-ix*4-iy*2-1)
428        end if
429 
430        neighbour = [ix,iy,iz]
431 
432        ! From indices_dense -> indices_coarse
433        curindices_coarse = curindices_dense(1:3) + neighbour(:)
434 
435        ! From indices_coarse -> ik_ibz in the coarse mesh
436        call get_kpt_from_indices_coarse(curindices_coarse,double_grid%maxcomp_coarse,&
437 &        double_grid%inttoik_coarse,double_grid%g0_coarse,double_grid%nbz_closedcoarse,ik_coarse,g0)
438 
439        ! Take into account a possible umklapp between k_dense and k_coarse
440        diffg0 = g0 - g01
441 
442        if (ANY(diffg0/=0)) then
443          ! WARNING works only with nspinor = 1 !!!
444          call calc_ceigr(diffg0,nfft,nspinor,Wfd_coarse%ngfft,ceigr)
445        end if
446 
447        do ib_dense = BSp%lomo_spin(spin), BSp%humo_spin(spin)
448          ! ur(ib_dense, ik_dense)
449          call wfd_sym_ur(Wfd_dense,Cryst,Kmesh_dense,ib_dense,ik_dense,spin,ur_dense)
450 
451          if (ANY(diffg0/=0)) then
452            !ur_kbz = ur_kbz*e(ig0r)
453            ur_dense(:) = ur_dense(:)*ceigr(:)
454          end if
455 
456          ! Uncomment for the complete overlap
457          !bstart = BSp%lomo_spin(spin); bstop = BSp%humo_spin(spin)
458 
459          ! Compute only dvv or dcc
460          if (ib_dense <= BSp%homo_spin(spin)) then
461            ! if ib_dense is a valence band => loop on valence bands
462            bstart = BSp%lomo_spin(spin); bstop = BSp%homo_spin(spin)
463          else
464            ! if ib_dense is a conduction band => loop on conduction bands
465            bstart = BSp%lumo_spin(spin); bstop = BSp%humo_spin(spin)
466          end if
467 
468          do ib_coarse = bstart, bstop
469            ! ur(ib_coarse, ik_coarse)
470            call wfd_sym_ur(Wfd_coarse,Cryst,Kmesh_coarse,ib_coarse,ik_coarse,spin,ur_coarse)
471 
472            ! ovlp = < u_{ib_coarse,ik_coarse} | u_{ib_dense,ik_dense} >
473            ovlp =  xdotc(nfft,ur_coarse,1,ur_dense,1)/nfft
474 
475            ! Filter too low values
476            if (ABS(ovlp) < threshold) ovlp = czero
477 
478            interpolator%overlaps(ib_coarse,ib_dense,ivertex,ik_dense,spin) = ovlp
479          end do ! ib_coarse
480 
481          !DBYG
482          ! write(std_out,*) "nb = ",neighbour
483          ! write(std_out,*) "(i1,i2,i3,j1,j2,j3) = ",curindices_dense
484          ! write(std_out,*) "ib = ",ib_dense
485          ! write(std_out,*) "Sum of dbb = ",REAL(SUM(GWPC_CONJG(interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin))*interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin))); call flush(std_out)
486          !ENDDBYG
487 
488        end do ! ib_dense
489      end do ! ivertex
490    end do ! ik_dense
491  end do ! spin
492 
493  ABI_FREE(ur_coarse)
494  ABI_FREE(ur_dense)
495  ABI_FREE(ceigr)
496 
497  ! Gather results on each node.
498  call xmpi_sum(interpolator%overlaps,Wfd_coarse%comm,ierr)
499 
500 end subroutine int_compute_overlaps

m_bseinterp/int_free_work [ Functions ]

[ Top ] [ m_bseinterp ] [ Functions ]

NAME

 int_free_work

FUNCTION

 Deallocate temporary arrays

INPUTS

OUTPUT

PARENTS

      m_hexc

CHILDREN

SOURCE

286 subroutine int_free_work(interpolator)
287 
288 
289 !This section has been created automatically by the script Abilint (TD).
290 !Do not modify the following lines by hand.
291 #undef ABI_FUNC
292 #define ABI_FUNC 'int_free_work'
293 !End of the abilint section
294 
295  implicit none
296 
297 !Arguments ---------------------------
298 !scalars
299  type(interpolator_t),intent(inout) :: interpolator
300 
301 !*****************************************************************************
302 
303  if( allocated(interpolator%btemp)) then
304    ABI_FREE(interpolator%btemp)
305  end if
306  if( allocated(interpolator%ctemp)) then
307    ABI_FREE(interpolator%ctemp)
308  end if
309 
310 end subroutine int_free_work

m_bseinterp/int_preprocess_tables [ Functions ]

[ Top ] [ m_bseinterp ] [ Functions ]

NAME

 int_preprocess_tables

FUNCTION

 Pre-process tables to improve interpolation technique

INPUTS

OUTPUT

PARENTS

      m_bseinterp

CHILDREN

SOURCE

523 subroutine int_preprocess_tables(interpolator,double_grid)
524 
525 
526 !This section has been created automatically by the script Abilint (TD).
527 !Do not modify the following lines by hand.
528 #undef ABI_FUNC
529 #define ABI_FUNC 'int_preprocess_tables'
530 !End of the abilint section
531 
532  implicit none
533 
534 !Argument ------------------------------------
535 !scalars
536  type(interpolator_t),intent(inout) :: interpolator
537  type(double_grid_t),intent(in) :: double_grid
538 !arrays
539 
540 !Local variables -----------------------------
541 !scalars
542  integer :: iorder,ik_dense,ik_coarse
543  integer :: ix,iy,iz,ineighbour,curdim, curj
544  real(dp) :: interp_factor
545 !arrays
546  integer :: allxyz(3),curindices_dense(6)
547  integer,allocatable :: curindex(:)
548 
549 !*********************************************
550  ABI_MALLOC(curindex,(double_grid%nbz_coarse))
551  curindex = 1
552 
553  interpolator%interp_factors = zero
554 
555  do ik_dense = 1,double_grid%nbz_dense
556 
557    ! From ik_ibz in the dense mesh -> indices_dense
558    iorder = double_grid%iktoint_dense(ik_dense)
559    !g01 = double_grid%g0_dense(:,iorder)
560 
561    ! From indices_dense -> indices_coarse
562    curindices_dense = double_grid%indices_dense(:,iorder)
563 
564    ik_coarse = double_grid%dense_to_coarse(ik_dense)
565 
566    ! Compute multi-linear interpolation factors
567    ! Loop over the neighbours
568    do ineighbour = 1,interpolator%nvert
569      !TODO helper function from [ix,iy,iz] -> ineighbour and vice versa
570      ix = (ineighbour-1)/4
571      iy = (ineighbour-ix*4-1)/2
572      iz = (ineighbour-ix*4-iy*2-1)
573      allxyz = [ix,iy,iz]
574      interp_factor = one
575      do curdim = 1,3
576        if (interpolator%method == BSE_INTERP_RL) then
577          cycle
578        else if(interpolator%method == BSE_INTERP_RL2) then
579          if (curdim /= 3) cycle
580        end if
581        curj = curindices_dense(3+curdim)
582        interp_factor = interp_factor*((allxyz(curdim)*(curj*1.0/double_grid%kmult(curdim)))&
583 &                               +((1-allxyz(curdim))*(1-(curj*1.0/double_grid%kmult(curdim)))))
584      end do
585      interpolator%interp_factors(ineighbour,curindex(ik_coarse)) = interp_factor
586    end do
587 
588    curindex(ik_coarse) = curindex(ik_coarse) + 1
589  end do
590 
591  ABI_FREE(curindex)
592 
593 end subroutine int_preprocess_tables

m_bseinterp/interpolator_free [ Functions ]

[ Top ] [ m_bseinterp ] [ Functions ]

NAME

 interpolator_free

FUNCTION

 Destroy the interpolator object in memory

INPUTS

OUTPUT

PARENTS

      m_hexc

CHILDREN

SOURCE

788 subroutine interpolator_free(interpolator)
789 
790 
791 !This section has been created automatically by the script Abilint (TD).
792 !Do not modify the following lines by hand.
793 #undef ABI_FUNC
794 #define ABI_FUNC 'interpolator_free'
795 !End of the abilint section
796 
797  implicit none
798 
799 !Arguments ---------------------------
800  type(interpolator_t),intent(inout) :: interpolator
801 
802 !*****************************************************************************
803 
804  if( allocated(interpolator%overlaps) ) then
805    ABI_FREE(interpolator%overlaps)
806  end if
807 
808  if (allocated(interpolator%corresp)) then
809    ABI_FREE(interpolator%corresp)
810  end if
811 
812  if (allocated(interpolator%interp_factors)) then
813    ABI_FREE(interpolator%interp_factors)
814  end if
815 
816  if( associated(interpolator%double_grid) ) then
817    nullify(interpolator%double_grid)
818  end if
819 
820 end subroutine interpolator_free

m_bseinterp/interpolator_init [ Functions ]

[ Top ] [ m_bseinterp ] [ Functions ]

NAME

 interpolator_init

FUNCTION

 Construct the interpolator object

INPUTS

OUTPUT

PARENTS

      m_hexc

CHILDREN

SOURCE

135 subroutine interpolator_init(interpolator, double_grid, Wfd_dense, Wfd_coarse, &
136 &    Kmesh_dense, Kmesh_coarse, BSp, Cryst, Psps, Pawtab, method)
137 
138 
139 !This section has been created automatically by the script Abilint (TD).
140 !Do not modify the following lines by hand.
141 #undef ABI_FUNC
142 #define ABI_FUNC 'interpolator_init'
143 !End of the abilint section
144 
145  implicit none
146 
147 !Arguments ---------------------------
148 !scalars
149  integer,intent(in) :: method
150  type(interpolator_t),intent(inout) :: interpolator
151  type(double_grid_t),intent(in),target :: double_grid
152  type(wfd_t),intent(inout) :: Wfd_dense, Wfd_coarse
153  type(kmesh_t),intent(in) :: Kmesh_dense, Kmesh_coarse
154  type(excparam),intent(in) :: BSp
155  type(crystal_t),intent(in) :: Cryst
156  type(pseudopotential_type),intent(in) :: Psps
157 !arrays
158  type(pawtab_type),intent(in) :: Pawtab(Cryst%ntypat*Wfd_coarse%usepaw)
159 
160 !Local variables ---------------------
161 !scalars
162  integer :: nsppol, nvert
163  integer :: maxnreh, nreh1, nreh2
164  integer :: mbandc, mbandd, nbzd
165  real(dp),parameter :: threshold = 0.1_dp
166 !arrays
167  character(len=500) :: msg
168 
169 !*****************************************************************************
170 
171  ABI_CHECK(Wfd_coarse%usepaw==0, "PAW not yet supported")
172  ABI_CHECK(BSp%nsppol==1, "nsppol != 1 not yet implemented")
173  ABI_CHECK(Wfd_coarse%nspinor==1, "nspinor != 1 not supported")
174 
175  ABI_UNUSED(Pawtab(1)%basis_size)
176  !paw_overlap(cprj1,cprj2,typat,pawtab,spinor_comm) result(onsite)
177 
178  interpolator%double_grid => double_grid
179  interpolator%mband_dense = Wfd_dense%mband
180  interpolator%mband_coarse = Wfd_coarse%mband
181  interpolator%method = method
182  interpolator%nsppol = BSp%nsppol
183 
184  SELECT CASE(method)
185  CASE (BSE_INTERP_YG)
186    nvert = 8
187  CASE (BSE_INTERP_RL2)
188    nvert = 2
189  CASE (BSE_INTERP_RL)
190    nvert = 1
191  CASE DEFAULT
192    write(msg,'(a,i0)') "Wrong interpolation method: ",method
193    MSG_ERROR(msg)
194  END SELECT
195 
196  interpolator%nvert = nvert
197 
198  mbandc = interpolator%mband_coarse
199  mbandd = interpolator%mband_dense
200  nbzd = double_grid%nbz_dense
201  nsppol = interpolator%nsppol
202  ABI_MALLOC(interpolator%overlaps,(mbandc,mbandd,nvert,nbzd,nsppol))
203 
204  call int_compute_overlaps(interpolator,double_grid, Wfd_dense, Wfd_coarse, Kmesh_dense, &
205 &   Kmesh_coarse, BSp, Cryst, Psps, Pawtab)
206 
207  ABI_MALLOC(interpolator%interp_factors,(nvert,double_grid%ndiv))
208 
209  call int_preprocess_tables(interpolator,double_grid)
210 
211  nreh1 = BSp%nreh(1)
212  nreh2 = nreh1; if(BSp%nsppol == 2) nreh2 = BSp%nreh(2)
213  maxnreh = MAX(nreh1,nreh2)
214 
215  ABI_MALLOC(interpolator%corresp,(maxnreh,interpolator%nvert,interpolator%nsppol))
216 
217  call int_compute_corresp(interpolator,BSp,double_grid)
218 
219 end subroutine interpolator_init

m_bseinterp/interpolator_normalize [ Functions ]

[ Top ] [ m_bseinterp ] [ Functions ]

NAME

 interpolator_normalize

FUNCTION

 Normalize the overlaps so that \sum_{ib} | d_{kk'}^{b,ib} | ^2 = 1

INPUTS

OUTPUT

PARENTS

      m_hexc

CHILDREN

SOURCE

720 subroutine interpolator_normalize(interpolator)
721 
722 
723 !This section has been created automatically by the script Abilint (TD).
724 !Do not modify the following lines by hand.
725 #undef ABI_FUNC
726 #define ABI_FUNC 'interpolator_normalize'
727 !End of the abilint section
728 
729  implicit none
730 
731 !Arguments ---------------------------
732 !scalars
733  type(interpolator_t),intent(inout) :: interpolator
734 
735 !Local variables ---------------------
736 !scalars
737  integer :: spin, ivertex, ib_dense, ik_dense
738  complex(gwpc) :: sum_ovlp
739 !arrays
740  complex(gwpc),allocatable :: overlaps(:)
741 
742 !*****************************************************************************
743 
744  ABI_MALLOC(overlaps,(interpolator%mband_coarse))
745  do spin = 1, interpolator%nsppol
746    do ivertex = 1, interpolator%nvert
747      do ib_dense = 1, interpolator%mband_dense
748        do ik_dense = 1, interpolator%double_grid%nbz_dense
749          overlaps(:) = interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin)
750          sum_ovlp = SQRT(REAL(SUM(GWPC_CONJG(overlaps(:))*overlaps(:))))
751          if (ABS(sum_ovlp) > tol6) then
752            overlaps(:) = overlaps(:)/sum_ovlp
753          else
754            overlaps(:) = czero
755          end if
756          interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin) = overlaps(:)
757        end do
758      end do
759    end do
760  end do
761  ABI_FREE(overlaps)
762 
763 end subroutine interpolator_normalize

m_haydock/int_compute_corresp [ Functions ]

[ Top ] [ m_haydock ] [ Functions ]

NAME

 int_compute_corresp

FUNCTION

INPUTS

 BSp<type(excparam)=The parameter for the Bethe-Salpeter run.
 grid <double_grid_t>=Correspondence between coarse and fine k-grid
 spin=Spin index.

OUTPUT

 corresp(Bsp%nreh(spin),8)= Correspondence between a transition on the
   coarse mesh and its i-th neighbour for i in [1,2,..,8].

 TODO:
  Some operations are faster if we allocate with shape (8,nreh(spin))

PARENTS

      m_bseinterp

CHILDREN

SOURCE

623 subroutine int_compute_corresp(interpolator,BSp,double_grid)
624 
625 
626 !This section has been created automatically by the script Abilint (TD).
627 !Do not modify the following lines by hand.
628 #undef ABI_FUNC
629 #define ABI_FUNC 'int_compute_corresp'
630 !End of the abilint section
631 
632  implicit none
633 
634 !Arguments ------------------------------------
635 !scalars
636  type(excparam),intent(in) :: BSp
637  type(double_grid_t),intent(in) :: double_grid
638  type(interpolator_t),intent(inout) :: interpolator
639 
640 !Local variables ------------------------------
641 !scalars
642  integer :: spin
643  integer :: itt,ik_dense,ik_coarse,iorder,it_coarse
644  integer :: ic,iv,ik_coarse0,it_coarse0,iovlp,ix,iy,iz
645 !arrays
646  integer :: curindices_dense(6),curindices_coarse(3),g0(3),g01(3),neighbour(3)
647 
648 !************************************************************************
649 
650  do spin=1,interpolator%nsppol
651    do itt=1,BSp%nreh_interp(spin)
652      ! From dense itt -> ik_dense, ic, iv
653      ik_dense = BSp%Trans_interp(itt,spin)%k
654      ic = BSp%Trans_interp(itt,spin)%c
655      iv = BSp%Trans_interp(itt,spin)%v
656 
657      ! From ik_dense -> indices_dense
658      iorder = double_grid%iktoint_dense(ik_dense)
659      g01 = double_grid%g0_dense(:,iorder)
660 
661      ! Index of the k-point in the coarse mesh.
662      ik_coarse0 = double_grid%dense_to_coarse(ik_dense)
663      it_coarse0 = BSp%vcks2t(iv,ic,ik_coarse0,spin)
664 
665      ! From indices_dense -> indices_coarse
666      curindices_dense = double_grid%indices_dense(:,iorder)
667 
668      ! Loop over the 8 neighbors.
669      do iovlp = 1,interpolator%nvert
670 
671        !TODO : helper function from [ix,iy,iz] -> iovlp and vice versa
672        if(interpolator%nvert > 1) then
673          ix = (iovlp-1)/4
674          iy = (iovlp-ix*4-1)/2
675          iz = (iovlp-ix*4-iy*2-1)
676        else
677          ix = (BSp%rl_nb-1)/4
678          iy = (BSp%rl_nb-ix*4-1)/2
679          iz = (BSp%rl_nb-ix*4-iy*2-1)
680        end if
681        neighbour = [ix,iy,iz]
682 
683        curindices_coarse = curindices_dense(1:3) + neighbour(:)
684 
685        ! From indices_coarse -> ik_ibz in the coarse mesh
686        call get_kpt_from_indices_coarse(curindices_coarse,double_grid%maxcomp_coarse,&
687 &        double_grid%inttoik_coarse,double_grid%g0_coarse,double_grid%nbz_closedcoarse,ik_coarse,g0)
688 
689        ! From ik_coarse, ic, iv to it_coarse
690        it_coarse = BSp%vcks2t(iv,ic,ik_coarse,spin)
691 
692        interpolator%corresp(it_coarse0,iovlp,spin) = it_coarse
693      end do
694    end do ! itt
695  end do
696 
697 end subroutine int_compute_corresp

m_haydock/interpolator_t [ Types ]

[ Top ] [ m_haydock ] [ Types ]

NAME

 interpolator_t

FUNCTION

  Store the overlap matrix elements needed for the interpolation of the BSE Hamiltonian

TODO

  Decide if we want to make the number of bands k-dependent.

SOURCE

 69  type,public :: interpolator_t
 70 
 71     integer :: nvert=8
 72     ! Number of vertices for interpolation
 73 
 74     integer :: method
 75     ! Interpolation method (YG or Rohlfing & Louie or ...)
 76 
 77     integer :: mband_dense, mband_coarse
 78     ! Max number of bands dense and coarse
 79 
 80     integer :: nsppol
 81     ! Number of spin channels
 82 
 83     integer, allocatable :: corresp(:,:,:)
 84     ! corresp(max_nreh,nvert,spin)
 85     ! it_coarse, idiv -> it_coarse (idiv-th neighbour)
 86 
 87     real(dp),allocatable :: interp_factors(:,:)
 88     ! interp_factors(nvert,ndiv)
 89     ! index_in_fine_box -> k-point in Trans_interp
 90 
 91     complex(gwpc),allocatable :: overlaps(:,:,:,:,:)
 92     ! Overlaps between dense and coarse mesh
 93     ! overlaps(mband_coarse,mband_dense,ivertex_coarse,double_grid%nkpt_dense,spin)
 94 
 95     complex(dpc),allocatable :: btemp(:), ctemp(:)
 96     ! Temporary arrays for work
 97 
 98     ! Pointers to datatypes that are already in memory
 99     type(double_grid_t),pointer :: double_grid => null()
100     ! Mapping between coarse and dense mesh
101 
102  end type interpolator_t
103 
104  public :: interpolator_init    ! Construct the object
105  public :: interpolator_free    ! Free memory
106  public :: interpolator_normalize ! Normalize the overlaps
107  public :: int_alloc_work       ! Alloc temp memory
108  public :: int_free_work        ! Free temp memory