TABLE OF CONTENTS


ABINIT/m_screen [ Modules ]

[ Top ] [ Modules ]

NAME

  m_screen

FUNCTION

   Screening object used in the BSE code.

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (MG)
  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

19 #if defined HAVE_CONFIG_H
20 #include "config.h"
21 #endif
22 
23 #include "abi_common.h"
24 
25 MODULE m_screen
26 
27  use defs_basis
28  use defs_abitypes
29  use m_xmpi
30  use m_hide_blas
31  use m_errors
32  use m_splines
33  use m_abicore
34  use m_kxc
35  use m_screening
36  use m_nctk
37  use m_sort
38 
39  use m_gwdefs,         only : GW_TOLQ0, czero_gw
40  use m_fstrings,       only : firstchar, endswith, strcat, itoa, sjoin
41  use m_numeric_tools,  only : print_arr
42  use m_geometry,       only : normv
43  use m_crystal,        only : crystal_t
44  use m_bz_mesh,        only : kmesh_t, get_BZ_item, has_bz_item
45  use m_gsphere,        only : gsphere_t
46  use m_vcoul,          only : vcoul_t
47  use m_io_screening,   only : hscr_free, hscr_io, read_screening, write_screening, hscr_print, &
48 &                             hscr_copy, hscr_t, hscr_bcast, hscr_from_file, ncname_from_id, em1_ncname, chi0_ncname
49  use m_ppmodel,        only : ppmodel_t, ppm_init, ppm_free, ppm_nullify, PPM_NONE, new_setup_ppmodel, ppm_symmetrizer
50 
51  implicit none
52 
53  private
54 
55  ! Flags defining the content of the %mat buffer in the fgg_t type.
56  integer,public,parameter :: MAT_NOTYPE         = 0
57  integer,public,parameter :: MAT_CHI0           = 1
58  integer,public,parameter :: MAT_CHI            = 2
59  integer,public,parameter :: MAT_EPSILON        = 3
60  integer,public,parameter :: MAT_INV_EPSILON    = 4
61  integer,public,parameter :: MAT_INV_EPSILON_M1 = 5
62  integer,public,parameter :: MAT_W              = 6
63  integer,public,parameter :: MAT_W_M1           = 7
64  !
65  ! Family vertex.
66  integer,public,parameter :: VTX_FAMILY_NONE  = 0   ! No vertex correction.
67  integer,public,parameter :: VTX_FAMILY_TDDFT = 1   ! TDDFT-based vertex.
68  integer,public,parameter :: VTX_FAMILY_ADA   = 2   ! ADA vertex.
69  !
70  ! Test charge or test particle.
71  integer,public,parameter :: VTX_TEST_CHARGE   = 0
72  integer,public,parameter :: VTX_TEST_PARTICLE = 1
73  !
74  ! Named constants for the frequency mesh.
75  integer,private,parameter :: WMESH_LINEAR    = 1
76  integer,private,parameter :: WMESH_GAUSS_LEG = 2
77  integer,private,parameter :: WMESH_TAN_GRID  = 3
78  !
79  ! Method used for the frequency integration.
80  integer,public,parameter ::  WINT_NONE    = 0
81  integer,public,parameter ::  WINT_PPMODEL = 1
82  integer,public,parameter ::  WINT_CONTOUR = 2
83  integer,public,parameter ::  WINT_AC      = 3
84  !
85  ! Parameters used for the model dielectric function.
86  integer,public,parameter :: MDL_NONE      = 0
87  integer,public,parameter :: MDL_BECHSTEDT = 1
88  !
89  ! Flags giving the status of the local buffers defined in fgg_t.
90  integer,private,parameter :: MAT_NODATA    = 0
91  integer,private,parameter :: MAT_ALLOCATED = 1
92  integer,private,parameter :: MAT_STORED    = 2
93  !
94  ! Flags giving the status of the local buffers defined in fgg_t.
95  integer,private,parameter :: FGG_QBZ_ISPOINTER  =1 ! Fgg_qbz is used to store the address in memory.
96  integer,private,parameter :: FGG_QBZ_ISALLOCATED=2 ! Fgg_qbz is used as an allocable array.

m_screen/em1_symmetrize_ip [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

  em1_symmetrize_ip

FUNCTION

  Symmetrizes the two-point function in G-space. Symmetrization is done
  inplace thorugh an auxiliary work array of dimension (npwc,npwc)

INPUTS

  nomega=All frequencies from 1 up to nomega are symmetrized.
  npwc=Number of G vectors in the symmetrized matrix.
  Gsph<gsphere_t>=data related to the G-sphere
  Qmesh<kmesh_t>=Structure defining the q-mesh used for Er.
  iq_bz=Index of the q-point in the BZ where epsilon^-1 is required.

SIDE EFFECTS

  epsm1(npwc,npwc,nomega)
   input:  filled with the matrix at the q-point that has to be symmetrized.
   output: symmetrised matrix.

NOTES

  In the present implementation we are not considering a possible umklapp vector G0 in the
  expression Sq = q+G0. Treating this case would require some changes in the G-sphere
  since we have to consider G-G0. The code however stops in sigma if a nonzero G0 is required
  to reconstruct the BZ.

  * Remember the symmetry properties of E
    If q_bz=Sq_ibz+G0:

    $ E_{SG1-G0,SG2-G0}(q_bz) = e^{+iS(G2-G1).\tau} E_{G1,G2)}(q)

    The invariance under exchange of the real space position E(1,2) = E(2,1) leads to:
    $ E_{-G2,-G1}(-q) = E_{G1,G2)

PARENTS

      m_screen

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

1566 subroutine em1_symmetrize_ip(iq_bz,npwc,nomega,Gsph,Qmesh,epsm1)
1567 
1568 
1569 !This section has been created automatically by the script Abilint (TD).
1570 !Do not modify the following lines by hand.
1571 #undef ABI_FUNC
1572 #define ABI_FUNC 'em1_symmetrize_ip'
1573 !End of the abilint section
1574 
1575  implicit none
1576 
1577 !Arguments ------------------------------------
1578 !scalars
1579  integer,intent(in) :: iq_bz,nomega,npwc
1580  type(gsphere_t),intent(in) :: Gsph
1581  type(kmesh_t),intent(in) :: Qmesh
1582 !arrays
1583  complex(gwpc),intent(inout) :: epsm1(npwc,npwc,nomega)
1584 
1585 !Local variables-------------------------------
1586 !scalars
1587  integer :: iw,g1,g2,isg1,isg2,iq_ibz,itim_q,isym_q,ierr
1588  logical :: q_isirred
1589  complex(gwpc) :: phmsg1t,phmsg2t_star
1590  character(len=500) :: msg
1591 !arrays
1592  real(dp) :: qbz(3)
1593  complex(gwpc),allocatable :: work(:,:)
1594 
1595 ! *********************************************************************
1596 
1597  ! * Get iq_ibz, and symmetries from iq_ibz.
1598  call get_BZ_item(Qmesh,iq_bz,qbz,iq_ibz,isym_q,itim_q,isirred=q_isirred)
1599 
1600  if (q_isirred) RETURN ! Nothing to do
1601 
1602  write(msg,'(a,f8.2,a)')" out of memory in work , requiring ",npwc**2*gwpc*b2Mb," Mb"
1603  ABI_STAT_MALLOC(work,(npwc,npwc), ierr)
1604  ABI_CHECK(ierr==0, msg)
1605 
1606 !$OMP PARALLEL DO PRIVATE(isg2,isg1,phmsg1t,phmsg2t_star,work) IF (nomega > 1)
1607  do iw=1,nomega
1608    do g2=1,npwc
1609      isg2 = Gsph%rottb(g2,itim_q,isym_q)
1610      phmsg2t_star = CONJG(Gsph%phmSGt(g2,isym_q))
1611      do g1=1,npwc
1612        isg1 = Gsph%rottb(g1,itim_q,isym_q)
1613        phmsg1t = Gsph%phmSGt(g1,isym_q)
1614        work(isg1,isg2) = epsm1(g1,g2,iw) * phmsg1t  * phmsg2t_star
1615      end do
1616    end do
1617    epsm1(:,:,iw) = work
1618  end do
1619 
1620  ABI_FREE(work)
1621  !
1622  ! Account for time-reversal ----------------------
1623  if (itim_q==2) then
1624 !$OMP PARALLEL DO IF (nomega > 1)
1625    do iw=1,nomega
1626      call sqmat_itranspose(npwc,epsm1(:,:,iw))
1627    end do
1628  end if
1629 
1630 end subroutine em1_symmetrize_ip

m_screen/em1_symmetrize_op [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

  em1_symmetrize_op

FUNCTION

  Symmetrizes the two-point function in G-space. Symmetrization is done outofplace.

INPUTS

  nomega=All frequencies from 1 up to nomega are symmetrized.
  npwc=Number of G vectors in the symmetrized matrix.
  Gsph<gsphere_t>=data related to the G-sphere
  Qmesh<kmesh_t>=Structure defining the q-mesh used for Er.
  iq_bz=Index of the q-point in the BZ where epsilon^-1 is required.
  in_epsm1(npwc,npwc,nomega)

OUTPUT

  out_epsm1(npwc,npwc,nomega)

NOTES

  In the present implementation we are not considering a possible umklapp vector G0 in the
  expression Sq = q+G0. Treating this case would require some changes in the G-sphere
  since we have to consider G-G0. The code however stops in sigma if a nonzero G0 is required
  to reconstruct the BZ.

  * Remember the symmetry properties of E
    If q_bz=Sq_ibz+G0:

    $ E_{SG1-G0,SG2-G0}(q_bz) = e^{+iS(G2-G1).\tau} E_{G1,G2)}(q)

    The invariance under exchange of the real space position E(1,2) = E(2,1) leads to:
    $ E_{-G2,-G1}(-q) = E_{G1,G2)

PARENTS

      m_screen

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

1675 subroutine em1_symmetrize_op(iq_bz,npwc,nomega,Gsph,Qmesh,in_epsm1,out_epsm1)
1676 
1677 
1678 !This section has been created automatically by the script Abilint (TD).
1679 !Do not modify the following lines by hand.
1680 #undef ABI_FUNC
1681 #define ABI_FUNC 'em1_symmetrize_op'
1682 !End of the abilint section
1683 
1684  implicit none
1685 
1686 !Arguments ------------------------------------
1687 !scalars
1688  integer,intent(in) :: iq_bz,nomega,npwc
1689  type(gsphere_t),target,intent(in) :: Gsph
1690  type(kmesh_t),intent(in) :: Qmesh
1691 !arrays
1692  complex(gwpc),intent(in) :: in_epsm1(npwc,npwc,nomega)
1693  complex(gwpc),intent(out) :: out_epsm1(npwc,npwc,nomega)
1694 
1695 !Local variables-------------------------------
1696 !scalars
1697  integer :: iw,g1,g2,isg1,isg2,iq_ibz,itim_q,isym_q
1698  logical :: q_isirred
1699  complex(gwpc) :: phmsg1t,phmsg2t_star
1700 !arrays
1701  real(dp) :: qbz(3)
1702 
1703 ! *********************************************************************
1704 
1705  ! * Get iq_ibz, and symmetries from iq_ibz.
1706  call get_BZ_item(Qmesh,iq_bz,qbz,iq_ibz,isym_q,itim_q,isirred=q_isirred)
1707 
1708  if (q_isirred) then
1709    out_epsm1 = in_epsm1; return
1710  end if
1711 
1712 ! grottb is a 1-1 mapping.
1713 !$OMP PARALLEL DO PRIVATE(isg1,isg2,phmsg1t,phmsg2t_star) COLLAPSE(2) IF (nomega > 1)
1714  do iw=1,nomega
1715    do g2=1,npwc
1716      isg2=Gsph%rottb(g2,itim_q,isym_q)
1717      phmsg2t_star = CONJG(Gsph%phmSGt(g2,isym_q))
1718      do g1=1,npwc
1719        isg1=Gsph%rottb(g1,itim_q,isym_q)
1720        phmsg1t = Gsph%phmSGt(g1,isym_q)
1721        out_epsm1(isg1,isg2,iw) = in_epsm1(g1,g2,iw) * phmsg1t * phmsg2t_star
1722      end do
1723    end do
1724  end do
1725  !
1726  ! Account for time-reversal ----------------------
1727  if (itim_q==2) then
1728 !$OMP PARALLEL DO IF (nomega > 1)
1729    do iw=1,nomega
1730      call sqmat_itranspose(npwc,out_epsm1(:,:,iw))
1731    end do
1732  end if
1733 
1734 end subroutine em1_symmetrize_op

m_screen/fgg_free_0D [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

 fgg_free_0D

FUNCTION

  Free memory.

PARENTS

      m_screen

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

427 subroutine fgg_free_0D(Fgg)
428 
429 
430 !This section has been created automatically by the script Abilint (TD).
431 !Do not modify the following lines by hand.
432 #undef ABI_FUNC
433 #define ABI_FUNC 'fgg_free_0D'
434 !End of the abilint section
435 
436  implicit none
437 
438 !Arguments ------------------------------------
439 !scalars
440  type(fgg_t),intent(inout) :: Fgg
441 ! *************************************************************************
442 
443  !@fgg_t
444  !complex
445  if (allocated(Fgg%mat))  then
446    ABI_FREE(Fgg%mat)
447  end if
448  Fgg%has_mat = MAT_NODATA
449 
450 end subroutine fgg_free_0D

m_screen/fgg_free_1D [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

 fgg_free_1D

FUNCTION

 Deallocate all the pointers in the structure that result to be associated.

 INPUT
  [keep_q(:)]=Optional logical mask used to select the q-points that are deallocated.

PARENTS

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

472 subroutine fgg_free_1D(Fgg,keep_q)
473 
474 
475 !This section has been created automatically by the script Abilint (TD).
476 !Do not modify the following lines by hand.
477 #undef ABI_FUNC
478 #define ABI_FUNC 'fgg_free_1D'
479 !End of the abilint section
480 
481  implicit none
482 
483 !Arguments ------------------------------------
484 !scalars
485  type(fgg_t),intent(inout) :: Fgg(:)
486  logical,optional,intent(in) :: keep_q(:)
487 
488 !Local variables ------------------------------
489 !scalars
490  integer :: iq_ibz
491  logical :: keep_it
492 ! *************************************************************************
493 
494  do iq_ibz=LBOUND(Fgg,DIM=1),UBOUND(Fgg,DIM=1)
495    keep_it = .FALSE.; if (PRESENT(keep_q)) keep_it = keep_q(iq_ibz)
496    if (.not. keep_it) call fgg_free_0D(Fgg(iq_ibz))
497  end do
498 
499 end subroutine fgg_free_1D

m_screen/fgg_init [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

 fgg_init

FUNCTION

 Initialize the structure allocating the memory and initializing the internal variables.

 INPUT
  npw
  nqlwl

SIDE EFFECTS

  Fgg, See function description

PARENTS

      m_screen

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

526 subroutine fgg_init(Fgg,npw,nomega,nqlwl)
527 
528 
529 !This section has been created automatically by the script Abilint (TD).
530 !Do not modify the following lines by hand.
531 #undef ABI_FUNC
532 #define ABI_FUNC 'fgg_init'
533 !End of the abilint section
534 
535  implicit none
536 
537 !Arguments ------------------------------------
538 !scalars
539  integer,intent(in) :: npw,nqlwl,nomega
540  type(fgg_t),intent(inout) :: Fgg
541 
542 !Local variables ------------------------------
543 !scalars
544  integer :: ierr
545 
546 ! *************************************************************************
547 
548  !@fgg_t
549  Fgg%nomega= nomega; Fgg%npw   = npw
550  Fgg%nqlwl = nqlwl
551 
552  if (npw>0.and.nomega>0) then
553    ABI_STAT_MALLOC(Fgg%mat,(npw,npw,nomega), ierr)
554    ABI_CHECK(ierr==0, "out of memory in Fgg%mat")
555    Fgg%has_mat = MAT_ALLOCATED
556  end if
557 
558 end subroutine fgg_init

m_screen/fgg_t [ Types ]

[ Top ] [ m_screen ] [ Types ]

NAME

 fgg_t

FUNCTION

  Object used to store F(G,G')(w) for a given q-point.

SOURCE

174  type,public :: fgg_t
175 
176   integer :: nomega
177   ! Number of frequencies.
178 
179   integer :: npw
180   ! Number of G vectors.
181 
182   integer :: nqlwl
183   ! Number of points for the treatment of the long wave-length limit.
184 
185   integer :: has_mat = MAT_NODATA
186   ! Flag giving the status of mat.
187 
188   !arrays
189   complex(gwpc),allocatable :: mat(:,:,:)
190   ! mat(npw,npw,nomega)
191   ! The component of the two-point function $F_{G,G',w}$ for a given q.
192 
193   !complex(dpc),allocatable :: head(:,:,:)
194   ! head(3,3,nomega)
195 
196   !complex(dpc),allocatable :: lwing(:,:,:)
197   ! lwing(3,npwe,nomega)
198   ! Lower wings
199 
200   !complex(dpc),allocatable :: uwing(:,:,:)
201   ! uwing(3,npwe,nomega)
202   ! Upper wings.
203 
204  end type fgg_t
205 
206  public :: fgg_init   ! Creation method.
207  public :: fgg_free   ! Free memory.

m_screen/screen_fgg_qbz_set [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

  screen_fgg_qbz_set

FUNCTION

  Helper function used to perform the setup W%Fgg_qbz setting also the internal
  flag that defines its status.

INPUTS

  W<screen_t>=Object describing the screened interaction.
  iq_bz=Index of the q-point in the BZ.
  nqlwl=Number of wings wanted.
  how= "Pointer" is a true pointer is wanted.
       "Allocated" if memory has to be allocated.

NOTES

  iq_bz and nqlwl are not used if how="Pointer".

PARENTS

      m_screen

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

589 subroutine screen_fgg_qbz_set(W,iq_bz,nqlwl,how)
590 
591 
592 !This section has been created automatically by the script Abilint (TD).
593 !Do not modify the following lines by hand.
594 #undef ABI_FUNC
595 #define ABI_FUNC 'screen_fgg_qbz_set'
596 !End of the abilint section
597 
598  implicit none
599 
600 !Arguments ------------------------------------
601 !scalars
602  integer,intent(in) :: iq_bz,nqlwl
603  character(len=*),intent(in) :: how
604  type(screen_t),intent(inout) :: W
605 
606 !Local variables ------------------------------
607 !scalars
608  !character(len=500) :: msg
609 
610 !************************************************************************
611 
612  W%fgg_qbz_idx = iq_bz ! Save the index of the q-point in the BZ.
613 
614  if (firstchar(how,(/"P"/)) ) then
615    ! We want a pointer.
616    select case (W%fgg_qbz_stat)
617    case (FGG_QBZ_ISALLOCATED)
618      call fgg_free_0D(W%Fgg_qbz)
619      ABI_DT_FREE(W%Fgg_qbz)
620      nullify(W%Fgg_qbz)
621      W%fgg_qbz_stat = FGG_QBZ_ISPOINTER
622 
623    case (FGG_QBZ_ISPOINTER)
624      ! Set it to null().
625      nullify(W%Fgg_qbz)
626 
627    case default
628      MSG_ERROR(sjoin("Wrong status:", itoa(W%fgg_qbz_stat)))
629    end select
630 
631  else if (firstchar(how,(/"A"/)) ) then
632    ! We want an allocatable array.
633 
634    select case (W%fgg_qbz_stat)
635    case (FGG_QBZ_ISPOINTER)
636      ! Allocate memory
637      nullify(W%Fgg_qbz)
638      ABI_DT_MALLOC(W%Fgg_qbz,)
639 
640      call fgg_init(W%Fgg_qbz,W%npw,W%nomega,nqlwl)
641      W%fgg_qbz_stat = FGG_QBZ_ISALLOCATED
642 
643    case (FGG_QBZ_ISALLOCATED)
644      W%Fgg_qbz%has_mat = MAT_ALLOCATED  ! STORED --> ALLOCATED
645 
646    case default
647      MSG_ERROR(sjoin("Wrong status:", itoa(W%fgg_qbz_stat)))
648    end select
649 
650  else
651    MSG_BUG(sjoin("Wrong how:", how))
652  end if
653 
654 end subroutine screen_fgg_qbz_set

m_screen/screen_free [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

 screen_free

FUNCTION

 Free the memory allocate in the datatype.

SIDE EFFECTS

 W<screen_t>=The data structure to be nullified.

OUTPUT

PARENTS

      bethe_salpeter

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

800 subroutine screen_free(W)
801 
802 
803 !This section has been created automatically by the script Abilint (TD).
804 !Do not modify the following lines by hand.
805 #undef ABI_FUNC
806 #define ABI_FUNC 'screen_free'
807 !End of the abilint section
808 
809  implicit none
810 
811 !Arguments ------------------------------------
812 !scalars
813  type(screen_t),intent(inout) :: W
814 
815 ! *************************************************************************
816 
817  !@screen_t
818  !
819  ! integer
820  if (allocated(W%gvec)) then
821    ABI_FREE(W%gvec)
822  end if
823  !
824  !real
825  if (allocated(W%ae_rhor)) then
826    ABI_FREE(W%ae_rhor)
827  end if
828 
829  if (allocated(W%qibz)) then
830    ABI_FREE(W%qibz)
831  end if
832 
833  if (allocated(W%qlwl)) then
834    ABI_FREE(W%qlwl)
835  end if
836  !
837  !complex
838  if (allocated(W%omega)) then
839    ABI_FREE(W%omega)
840  end if
841 
842  ! logical
843  if (allocated(W%keep_q)) then
844    ABI_FREE(W%keep_q)
845  end if
846  !
847  ! types ------------------------------------------
848  !
849  ! Here be careful with dangling pointers.
850  ! First Fgg_qbz that might point to one of the %Fgg then %Fgg.
851  select case (W%fgg_qbz_stat)
852  case (FGG_QBZ_ISALLOCATED)
853    call fgg_free_0D(W%Fgg_qbz)
854    ABI_DT_FREE(W%Fgg_qbz)
855    nullify(W%Fgg_qbz)
856    W%fgg_qbz_stat=FGG_QBZ_ISPOINTER
857 
858  case (FGG_QBZ_ISPOINTER)
859    nullify(W%Fgg_qbz)
860    W%fgg_qbz_stat=FGG_QBZ_ISPOINTER
861 
862  case default
863    continue
864  end select
865  !
866  ! Free the Fgg matrices.
867  if (associated(W%Fgg)) then
868    call fgg_free(W%Fgg)
869    ABI_DT_FREE(W%Fgg)
870  end if
871  !
872  ! Free the plasmon pole tables.
873  call ppm_free(W%PPm)
874 
875 end subroutine screen_free

m_screen/screen_ihave_fgg [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

  screen_ihave_fgg

FUNCTION

  Inquire the processor whether it has a particular F_{GG')(q) and with which status.

INPUTS

  W<screen_t>=Object describing the screened interaction.
  iq_ibz=k-point index
  [how]=string defining which status is checked. By default the function returns
     .TRUE. if the wave is either MAT_ALLOCATED or MAT_STORED.
     Possible mutually exclusive values: "Allocated", "Stored".
     Only the first character is checked (no case-sensitive)

NOTES

   A zero index can be used to inquire the status of the full set of q-points.

PARENTS

CHILDREN

SOURCE

683 pure function screen_ihave_fgg(W,iq_ibz,how)
684 
685 
686 !This section has been created automatically by the script Abilint (TD).
687 !Do not modify the following lines by hand.
688 #undef ABI_FUNC
689 #define ABI_FUNC 'screen_ihave_fgg'
690 !End of the abilint section
691 
692  implicit none
693 
694 !Arguments ------------------------------------
695 !scalars
696  integer,intent(in) :: iq_ibz
697  logical :: screen_ihave_fgg
698  character(len=*),optional,intent(in) :: how
699  type(screen_t),intent(in) :: W
700 
701 !Local variables ------------------------------
702 !scalars
703  integer :: ii
704 !arrays
705  integer :: check(2)
706 
707 !************************************************************************
708 
709  check = (/MAT_ALLOCATED, MAT_STORED/)
710  if (PRESENT(how)) then
711    if (firstchar(how,(/"A","a"/))) check = (/MAT_ALLOCATED, MAT_ALLOCATED/)
712    if (firstchar(how,(/"S","s"/))) check = (/MAT_STORED, MAT_STORED/)
713  end if
714 
715  if (iq_ibz>0) then
716    screen_ihave_fgg = (W%Fgg(iq_ibz)%has_mat==check(1) .or.&
717                        W%Fgg(iq_ibz)%has_mat==check(2) )
718  else
719    ! check the status of the full set of q-tables.
720    screen_ihave_fgg=.TRUE.
721    do ii=1,W%nqibz
722      screen_ihave_fgg = screen_ihave_fgg .and. &
723                      (W%Fgg(ii)%has_mat==check(1) .or.&
724                       W%Fgg(ii)%has_mat==check(2) )
725    end do
726  end if
727 
728 end function screen_ihave_fgg

m_screen/screen_info_print [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

  screen_info_print

FUNCTION

  Printout of screen_info_t.

INPUTS

  W_info<screen_info_t>=The structure.
  [unit]=Unit number for output
  [prtvol]=Verbosity level
  [mode_paral]=Either "COLL" or "PERS"
  [header]=String to be printed as header for additional info.

OUTPUT

  Only printing

PARENTS

      m_screen

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

346 subroutine screen_info_print(W_info,header,unit,mode_paral,prtvol)
347 
348 
349 !This section has been created automatically by the script Abilint (TD).
350 !Do not modify the following lines by hand.
351 #undef ABI_FUNC
352 #define ABI_FUNC 'screen_info_print'
353 !End of the abilint section
354 
355  implicit none
356 
357 !Arguments ------------------------------------
358 !scalars
359  integer,optional,intent(in) :: unit,prtvol
360  character(len=4),optional,intent(in) :: mode_paral
361  character(len=*),optional,intent(in) :: header
362  type(screen_info_t),intent(in) :: W_info
363 
364 !Local variables-------------------------------
365  integer :: my_unt,my_prtvol
366  character(len=4) :: my_mode
367  character(len=500) :: msg
368 ! *********************************************************************
369 
370  !@screen_info_t
371  my_unt   =std_out; if (PRESENT(unit      )) my_unt   =unit
372  my_prtvol=0      ; if (PRESENT(prtvol    )) my_prtvol=prtvol
373  my_mode  ='COLL' ; if (PRESENT(mode_paral)) my_mode  =mode_paral
374 
375  msg=' ==== Info on the screen_info_t% object ==== '
376  if (PRESENT(header)) msg=' ==== '//TRIM(ADJUSTL(header))//' ==== '
377  call wrtout(my_unt,msg,my_mode)
378 
379 !integer
380  write(msg,'(a,i3)')" mat_type    ",W_info%mat_type
381  call wrtout(my_unt,msg,my_mode)
382  write(msg,'(a,i3)')" vtx_family  ",W_info%vtx_family
383  call wrtout(my_unt,msg,my_mode)
384  write(msg,'(a,i3)')" invalid_freq",W_info%invalid_freq
385  call wrtout(my_unt,msg,my_mode)
386  write(msg,'(a,i3)')" ixc         ",W_info%ixc
387  call wrtout(my_unt,msg,my_mode)
388  write(msg,'(a,i3)')" use_ada     ",W_info%use_ada
389  call wrtout(my_unt,msg,my_mode)
390  write(msg,'(a,i3)')" use_mdf     ",W_info%use_mdf
391  call wrtout(my_unt,msg,my_mode)
392  write(msg,'(a,i3)')" use_ppm     ",W_info%use_ppm
393  call wrtout(my_unt,msg,my_mode)
394  write(msg,'(a,i3)')" vtx_test    ",W_info%vtx_test
395  call wrtout(my_unt,msg,my_mode)
396  write(msg,'(a,i3)')" wint_method ",W_info%wint_method
397  call wrtout(my_unt,msg,my_mode)
398 
399 !real
400  write(msg,'(a,f8.3)')" ada_kappa   ",W_info%ada_kappa
401  call wrtout(my_unt,msg,my_mode)
402  write(msg,'(a,f8.3)')" eps_inf     ",W_info%eps_inf
403  call wrtout(my_unt,msg,my_mode)
404  write(msg,'(a,f8.3)')" drude_plsmf ",W_info%drude_plsmf
405  call wrtout(my_unt,msg,my_mode)
406 
407 end subroutine screen_info_print

m_screen/screen_info_t [ Types ]

[ Top ] [ m_screen ] [ Types ]

NAME

 screen_info_t

FUNCTION

  Container storing the parameters used to initialize a screen_t datatype or to
  calculate a new SCR file from the SUSC file containing the independent-particle
  polarizability.

NOTES

  The list of parameters passed to screening_init is copied in W%Info.
  At present there is no need to provide a copy method since the structure does
  not contain pointers but such a method must be defined and used if the
  dynamic entities are added to the datatype.

SOURCE

118 type,public :: screen_info_t
119 
120   integer :: mat_type = MAT_NOTYPE
121   ! Matrix identifier. See MAT_* flags.
122 
123   integer :: vtx_family = VTX_FAMILY_NONE
124   ! Vertex correction family.
125 
126   integer :: invalid_freq=0
127   ! Sets the procedure to follow when a PPM frequency is invalid (negative or imaginary), see input variable gw_invalid_freq
128 
129   integer :: ixc=0
130   ! XC functional used for the TDDFT-based vertex.
131 
132   integer :: use_ada=0
133   ! >0 if ADA vertex is used.
134 
135   integer :: use_mdf = MDL_NONE
136   ! >0 if model dielectric function is used.
137 
138   integer :: use_ppm = PPM_NONE
139   ! >0 if ppmodel is used.
140 
141   integer :: vtx_test = VTX_TEST_CHARGE
142   ! test charge or test particle.
143 
144   integer :: wint_method = WINT_NONE
145   ! Defines the frequency integration technique. See WIN_ flags.
146   ! NOTE that this flag can be changed at run time. For example
147   ! one can switch from the CD to the PPm if the ppmodel parameters are in memory
148 
149   real(dp) :: ada_kappa = 2.1_dp
150   ! Inverse smearing length used for ADA.
151 
152   real(dp) :: eps_inf = 12.0_dp
153   ! Dielectric constant used for the model dielectric function.
154 
155   real(dp) :: drude_plsmf = zero
156   ! Drude plasma frequency used for PPmodel 1.
157 
158 end type screen_info_t
159 
160  public :: screen_info_print

m_screen/screen_init [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

  screen_init

FUNCTION

  Initialize basic dimensions and the important (small) arrays in an Epsilonm1_results data type
  starting from a file containing either epsilon^{-1} (_SCR) or chi0 (_SUSC).

INPUTS

  W_Info<screen_info_t>=The list of parameters used to construct the screen function.
  Cryst<crystal_t>=Info on the unit cell.
  Qmesh<kmesh_t>=Info on the Q-mesh.
  Gsph<gsphere_t>=Info on the plane-wave basis set used for the two-point function.
  Vcp<vcoul_t>=Structure gathering data on the Coulombian interaction.
  ifname=The name of the external file used to read the matrix.
  id_required=Identifier used to specify the type of two-point function that is wanted.
  iomode=Option defining the file format of the external file.
  mqmem=0 for out-of-core solution, /=0 if entire matrix has to be stored in memory.
  npw_asked=Number of G-vector to be used in the calculation, if <=0 use Max allowed number.
  ngfftf(18)=Info on the (fine) mesh used for the density.
  nfftf_tot=Total number of point in the FFT mesh for ae_rhor
  nsppol=Numer of independent spin polarizations.
  nspden=Number of spin density components in ae_rhor
  ae_rhor(nfftf_tot,nspden)
  prtvol=Verbosity level.
  comm=MPI communicator.

OUTPUT

  W<screen_t>=The structure initialized with basic dimensions and arrays.

PARENTS

      bethe_salpeter

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

 918 subroutine screen_init(W,W_Info,Cryst,Qmesh,Gsph,Vcp,ifname,mqmem,npw_asked,&
 919 &  iomode,ngfftf,nfftf_tot,nsppol,nspden,ae_rhor,prtvol,comm)
 920 
 921 
 922 !This section has been created automatically by the script Abilint (TD).
 923 !Do not modify the following lines by hand.
 924 #undef ABI_FUNC
 925 #define ABI_FUNC 'screen_init'
 926 !End of the abilint section
 927 
 928  implicit none
 929 
 930 !Arguments ------------------------------------
 931 !scalars
 932  integer,intent(in) :: mqmem,iomode,npw_asked,comm,prtvol,nsppol
 933  integer,intent(in) :: nfftf_tot,nspden
 934  character(len=fnlen),intent(in) :: ifname
 935  type(crystal_t),intent(in) :: Cryst
 936  type(gsphere_t),intent(in) :: Gsph
 937  type(vcoul_t),intent(in) :: Vcp
 938  type(kmesh_t),intent(in) :: Qmesh
 939  type(screen_info_t),intent(in) :: W_Info
 940  type(screen_t),intent(out) :: W
 941 !arrays
 942  integer,intent(in) :: ngfftf(18)
 943  real(dp),intent(in) :: ae_rhor(nfftf_tot,nspden)
 944 
 945 !Local variables-------------------------------
 946 !scalars
 947  integer,parameter :: master=0
 948  integer :: option_test,approx_type,ixc_required,id_required !,nkxc
 949  integer :: fform,my_rank,mat_type_read
 950  integer :: nqibz,nomega,iq_ibz,npw,nqlwl
 951  integer :: nI,nJ,iq_bz,mdf_type,ppmodel,ierr
 952  integer :: iw,qsort,ii
 953  real(dp) :: eps_inf,drude_plsmf
 954  logical :: deallocate_Fgg,found,from_file,is_qeq0,remove_dgg !only_one_kpt,
 955  character(len=500) :: msg
 956  character(len=fnlen) :: sus_fname,scr_fname
 957  character(len=nctk_slen) :: varname
 958  type(hscr_t) :: Hscr
 959 !arrays
 960  integer :: g0(3),iperm(Qmesh%nibz)
 961  real(dp) :: wt_list(Qmesh%nibz)
 962  !complex(gwpc),ABI_CONTIGUOUS pointer :: em1_ggw(:,:,:)
 963 
 964 ! *********************************************************************
 965 
 966  DBG_ENTER("COLL")
 967 
 968  ABI_UNUSED(nsppol)
 969 
 970  !@screen_t
 971  my_rank = xmpi_comm_rank(comm)
 972 
 973  call screen_nullify(W)
 974  !
 975  ! Initialize basic parameters ----------------------------------------
 976  !
 977  !@screen_info_t
 978  W%Info = W_info ! Copy basic parameters.
 979  call screen_info_print(W%Info,header="W info",unit=std_out)
 980 
 981  id_required  = W_Info%mat_type
 982  approx_type  = W_Info%vtx_family
 983  option_test  = W_Info%vtx_test
 984  ixc_required = W_Info%ixc
 985  varname = ncname_from_id(id_required)
 986 
 987  !if (ALL(id_required /= (/MAT_W_M1, MAT_W/)) ) then
 988  if (all(id_required /= [MAT_INV_EPSILON])) then
 989    MSG_ERROR(sjoin("id_required:",itoa(id_required),"not available"))
 990  end if
 991 
 992  ! This part must be rationalized.
 993  remove_dgg = (id_required==MAT_W_M1)
 994 
 995  if (W%Info%use_mdf == MDL_NONE) W%fname = ifname
 996  W%nI = 1; W%nJ = 1
 997  !
 998  ! The Q-point sampling is inited from Qmesh.
 999  W%nqibz=Qmesh%nibz
1000  ABI_MALLOC(W%qibz,(3,W%nqibz))
1001  W%qibz= Qmesh%ibz
1002 
1003  W%mqmem=mqmem!;  W%mqmem = 0
1004  if (W%mqmem/=0) W%mqmem=W%nqibz
1005 
1006  ABI_MALLOC(W%keep_q,(W%nqibz))
1007  W%keep_q=.TRUE.; if (W%mqmem == 0) W%keep_q = .False.
1008 
1009  if (W%mqmem/=0 .and. W%mqmem<W%nqibz) then
1010    ! Keep in memory the most representative q-points.
1011    W%keep_q=.FALSE.
1012    wt_list = Qmesh%wt; iperm = (/(ii,ii=1,Qmesh%nibz)/)
1013    call sort_dp(Qmesh%nibz,wt_list,iperm,tol12)
1014    do qsort=Qmesh%nibz,Qmesh%nibz-mqmem+1,1
1015      iq_ibz = iperm(qsort)
1016      W%keep_q(iq_ibz) = .TRUE.
1017    end do
1018  end if
1019  W%fgg_qbz_idx = 0
1020 
1021  W%iomode = iomode
1022  W%prtvol = prtvol
1023 
1024  W%has_ppmodel=0; if (W%Info%use_ppm/=PPM_NONE) W%has_ppmodel=1
1025  !
1026  ! Copy the AE density for the model dielectric function or for the vertex corrections.
1027  W%nspden     = nspden
1028  W%ngfftf     = ngfftf
1029  W%nfftf_tot  = nfftf_tot
1030 
1031  ABI_MALLOC(W%ae_rhor,(nfftf_tot,nspden))
1032  W%ae_rhor = ae_rhor
1033 
1034  deallocate_Fgg=.FALSE.
1035 
1036  W%has_fgg=0; if ( ANY(W%Info%wint_method == (/WINT_CONTOUR, WINT_AC/)) ) W%has_fgg=1
1037 
1038  if (W%has_fgg>0 .and. W%has_ppmodel>0) then
1039    MSG_WARNING("Both PPmodel tables and F_(GG')(q,w) are stored in memory")
1040  end if
1041 
1042  ! G-sphere for W and Sigma_c is initialized from the SCR file.
1043  !%% only_one_kpt=.FALSE.
1044  !%% call gsph_init(Gsph,Cryst,Sigp%npwc,gvec=Er%gvec)
1045  !%% deallocate(gvec_p)
1046 
1047  ! Default values used if external file is not read.
1048  nqlwl  = 0
1049  nomega = 1
1050  W%npw  = npw_asked
1051 
1052  ! Model dielectric function does not require any external file.
1053  from_file = (W%Info%use_mdf==MDL_NONE)
1054 
1055  if (from_file) then
1056    ! Open file and check its content.
1057    if (endswith(W%fname, ".nc")) W%iomode = IO_MODE_ETSF
1058 
1059    call hscr_from_file(hscr,W%fname,fform,comm)
1060    ! Echo of the header
1061    if (my_rank == master .and. W%prtvol>0) call hscr_print(hscr)
1062 
1063    ! Communicate the header and copy basic parameters.
1064    !call hscr_bcast(Hscr,master,my_rank,comm)
1065    mat_type_read = Hscr%id
1066    nqlwl         = Hscr%nqlwl
1067    nomega        = Hscr%nomega
1068  end if
1069 
1070  W%nqlwl  = nqlwl
1071  W%nomega = nomega
1072 
1073  ABI_MALLOC(W%qlwl,(3,W%nqlwl))
1074  ABI_MALLOC(W%omega,(W%nomega))
1075 
1076  if (from_file) then
1077    W%qlwl  = Hscr%qlwl
1078    W%omega = Hscr%omega
1079    !
1080    ! G-vectors.
1081    W%npw=Hscr%npwe
1082    if (npw_asked>0) then
1083      if (npw_asked>Hscr%npwe) then
1084        write(msg,'(a,i8,2a,i8)')&
1085         'The number of G-vectors saved on file is less than the value required = ',npw_asked,ch10,&
1086         'Calculation will proceed with the Max available npw = ',Hscr%npwe
1087        MSG_WARNING(msg)
1088      else
1089        W%npw=npw_asked ! Redefine the no. of G"s for W.
1090        write(msg,'(a,i8,2a,i8)')&
1091         'The Number of G-vectors saved on file is larger than the value required = ',npw_asked,ch10,&
1092         'Calculation will proceed with npw = ',W%npw
1093        MSG_COMMENT(msg)
1094      end if
1095    end if
1096    !
1097    ! consistency check on G-vectors and q-points.
1098    if ( ANY(Hscr%gvec(:,1:W%npw) /= Gsph%gvec(:,1:W%npw)) ) then
1099      !write(std_out) W%gvec, Gsph%gvec
1100      MSG_ERROR("Hscr%gvec /= Gsph%gvec(1:W%npw)")
1101    end if
1102    ABI_CHECK(Hscr%nqibz==Qmesh%nibz,"Mismatch in the number of q-points")
1103    ierr=0
1104    do iq_ibz=1,Hscr%nqibz
1105      if (ANY( ABS(Qmesh%ibz(:,iq_ibz) - Hscr%qibz(:,iq_ibz)) > tol6) ) then
1106        ierr=ierr+1
1107        write(std_out,'(i0,2(3f7.3,1x))')iq_ibz,Qmesh%ibz(:,iq_ibz),Hscr%qibz(:,iq_ibz)
1108      end if
1109    end do
1110    ABI_CHECK(ierr==0,"Wrong ordering in q-point list, aborting now")
1111  end if
1112 
1113  ABI_MALLOC(W%gvec,(3,W%npw))
1114  W%gvec = Gsph%gvec(:,1:W%npw)
1115 
1116  ! Frequency mesh.
1117  W%nomega_r=1; W%nomega_i=0
1118  if (W%nomega==2) then
1119    W%nomega_r=1; W%nomega_i=1
1120  else
1121    ! Real frequencies are packed in the first locations.
1122    W%nomega_r=1
1123    do iw=1,W%nomega
1124      if (DBLE(W%omega(iw))>0.001*Ha_eV) W%nomega_r=iw
1125    end do
1126    W%nomega_i=W%nomega-W%nomega_r
1127  end if
1128  !
1129  ! ------------------------------ Initialization completed --------------------------------
1130  !
1131  ! Just to keep the code readable.
1132  npw    = W%npw
1133  nqibz  = W%nqibz
1134  nomega = W%nomega
1135  nI     = W%ni
1136  nJ     = W%nj
1137  ABI_DT_MALLOC(W%Fgg,(nqibz))
1138 
1139  if (from_file) then
1140 
1141    ! Read the ab-initio em1 from file.
1142    select case (mat_type_read)
1143    case (MAT_INV_EPSILON)
1144      call wrtout(std_out,strcat("Em1 will be initialized from SCR file: ",W%fname),"COLL")
1145      !
1146    case  (MAT_CHI0)  ! Write new SCR file.
1147      MSG_ERROR("Not coded yet")
1148      sus_fname = W%fname; scr_fname="TESTING_SUS2SCR"
1149 
1150      W%fname = scr_fname  ! Change the name of the file associated to W.
1151      !
1152    case default
1153      write(msg,'(a,i0)')" Unsupported conversion from mat_type ",mat_type_read
1154      MSG_ERROR(msg)
1155    end select
1156 
1157    do iq_ibz=1,nqibz
1158      if (.not.W%keep_q(iq_ibz)) then
1159         call wrtout(std_out,strcat("Skipping iq_ibz: ",itoa(iq_ibz)),"COLL")
1160         CYCLE
1161      end if
1162      !
1163      nqlwl=0; is_qeq0=(normv(W%qibz(:,iq_ibz),Cryst%gmet,'G')<GW_TOLQ0)
1164      if (is_qeq0) nqlwl=W%nqlwl
1165      !
1166      ! Allocate F_{GG'}(w)
1167      call fgg_init(W%Fgg(iq_ibz),npw,nomega,nqlwl)
1168      !
1169      ! Read data from file (use MPI-IO if possible)
1170      if (W%iomode /= IO_MODE_ETSF .and. xmpi_mpiio==1) then
1171        call wrtout(std_out,ABI_FUNC//"read_screening with MPI_IO","COLL")
1172        call read_screening(varname,W%fname,npw,1,nomega,W%Fgg(iq_ibz)%mat,IO_MODE_MPI,comm,iqiA=iq_ibz)
1173      else
1174        call read_screening(varname,W%fname,npw,1,nomega,W%Fgg(iq_ibz)%mat,W%iomode,comm,iqiA=iq_ibz)
1175      end if
1176      W%Fgg(iq_ibz)%has_mat = MAT_STORED
1177      ! W contains Em1 and is ready to use.
1178    end do
1179 
1180  else
1181    !
1182    ! Model dielectric function. Only epsm-1 is supported here.
1183    call wrtout(std_out," Calculating model dielectric function... ","COLL")
1184    ABI_CHECK(W%nomega==1,"Cannot use nomega>1 in model dielectric function")
1185    !
1186    do iq_ibz=1,nqibz
1187      !
1188      if (.not.W%keep_q(iq_ibz)) CYCLE
1189      !
1190      ! The wings are not used here.
1191      nqlwl=0; is_qeq0= (normv(W%qibz(:,iq_ibz),Cryst%gmet,'G')<GW_TOLQ0)
1192      !
1193      ! Calculate the model. Note that mdielf awaits an index in the BZ.
1194      found = has_bz_item(Qmesh,Qmesh%ibz(:,iq_ibz),iq_bz,g0)
1195      if (.not.found.or.ANY(g0/=0)) then
1196        MSG_ERROR("Problem in retrieving ibz point")
1197      end if
1198      !
1199      ! Allocate F_{GG'}(w).
1200      call fgg_init(W%Fgg(iq_ibz),npw,nomega,nqlwl)
1201 
1202      eps_inf  =  W%Info%eps_inf
1203      mdf_type =  W%Info%use_mdf
1204      !em1_ggw  => W%Fgg(iq_ibz)%mat
1205 
1206      ! Construct W TODO check the new implementation.
1207      call screen_mdielf(iq_bz,npw,nomega,mdf_type,eps_inf,Cryst,Qmesh,Vcp,Gsph,&
1208 &      nspden,nfftf_tot,ngfftf,ae_rhor,"EM1",W%Fgg(iq_ibz)%mat,comm)
1209 
1210      W%Fgg(iq_ibz)%has_mat = MAT_STORED
1211 
1212      if (W%prtvol > 0) then
1213        do iw=1,nomega
1214          write(msg,'(a,i3,a,i4,a)')'  Model symmetrical e^{-1} (q=',iq_ibz,', omega=',iw,', G,G'')'
1215          call wrtout(std_out,msg,'COLL')
1216          call print_arr(W%Fgg(iq_ibz)%mat(:,:,iw))
1217        end do
1218      end if
1219      !
1220    end do ! iq_ibz
1221    !
1222  end if
1223  !
1224  ! Init plasmon-pole parameters.
1225  if (W%has_ppmodel>0) then
1226    MSG_WARNING("Calculating PPmodel parameters")
1227    ppmodel = W%Info%use_ppm; drude_plsmf = W%Info%drude_plsmf
1228    call ppm_init(W%PPm,W%mqmem,W%nqibz,W%npw,ppmodel,drude_plsmf,W%Info%invalid_freq)
1229 
1230    do iq_ibz=1,nqibz
1231      if (screen_ihave_fgg(W,iq_ibz,how="Stored")) then
1232        !em1_ggw => W%Fgg(iq_ibz)%mat
1233        call new_setup_ppmodel(W%PPm,iq_ibz,Cryst,Qmesh,npw,nomega,W%omega,W%Fgg(iq_ibz)%mat,&
1234           nfftf_tot,Gsph%gvec,ngfftf,W%ae_rhor(:,1))
1235      end if
1236    end do
1237  end if
1238  !
1239  ! Deallocate Fgg if the matrices are not needed anymore.
1240  if (deallocate_Fgg) then
1241    call screen_fgg_qbz_set(W,0,0,"Pointer") ! Avoid dangling pointer.
1242    call fgg_free(W%Fgg,keep_q=W%keep_q)
1243  end if
1244 
1245  if (from_file) call hscr_free(Hscr)
1246 
1247  DBG_EXIT("COLL")
1248 
1249 end subroutine screen_init

m_screen/screen_nullify [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

 screen_nullify

FUNCTION

 Initialize the pointers to null()

SIDE EFFECTS

 W<screen_t>=The data structure to be nullified.

PARENTS

      bethe_salpeter,m_screen

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

751 subroutine screen_nullify(W)
752 
753 
754 !This section has been created automatically by the script Abilint (TD).
755 !Do not modify the following lines by hand.
756 #undef ABI_FUNC
757 #define ABI_FUNC 'screen_nullify'
758 !End of the abilint section
759 
760  implicit none
761 
762 !Arguments ------------------------------------
763 !scalars
764  type(screen_t),intent(inout) :: W
765 ! *************************************************************************
766 
767  !@screen_t
768  !
769  !types
770  nullify(W%Fgg_qbz); W%fgg_qbz_stat=FGG_QBZ_ISPOINTER ! Nedded since the initial status is undefined.
771  nullify(W%Fgg)
772 
773  call ppm_nullify(W%PPm)
774 
775 end subroutine screen_nullify

m_screen/screen_symmetrizer [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

  screen_symmetrizer

FUNCTION

  Modify the status of the object so that the symmetrized component F(q_bz)_GG' is calculated
  (if needed) and is made available in the internal buffer. This routine must be called before
  performing any operation that involves the symmetrized component of the two-point function.

INPUTS

  W<screen_t>=Data structure used to represent the two-point function
  Cryst<crystal_t>=Info on the crystal structure.
  Gsph<gsphere_t>=data related to the G-sphere
  Qmesh<kmesh_t>=Structure defining the q-mesh used for sample the BZ.
  iq_bz=Index of the q-point in the BZ where F(q_bz)_GG' is wanted.

OUTPUT

SIDE EFFECTS

   W%PPm
   W%Fgg_qbz

NOTES

  In the present implementation we are not considering a possible umklapp vector G0 in the
  expression Sq = q+G0. Treating this case would require some changes in the G-sphere
  since we have to consider G-G0. The code however stops in sigma if a nonzero G0 is required
  to reconstruct the BZ.

PARENTS

      exc_build_block

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

1290 subroutine screen_symmetrizer(W,iq_bz,Cryst,Gsph,Qmesh,Vcp)
1291 
1292 
1293 !This section has been created automatically by the script Abilint (TD).
1294 !Do not modify the following lines by hand.
1295 #undef ABI_FUNC
1296 #define ABI_FUNC 'screen_symmetrizer'
1297 !End of the abilint section
1298 
1299  implicit none
1300 
1301 !Arguments ------------------------------------
1302 !scalars
1303  integer,intent(in) :: iq_bz
1304  type(screen_t),intent(inout) :: W
1305  type(crystal_t),intent(in) :: Cryst
1306  type(gsphere_t),intent(in) :: Gsph
1307  type(kmesh_t),intent(in) :: Qmesh
1308  type(vcoul_t),intent(in) :: Vcp
1309 
1310 !Local variables-------------------------------
1311 !scalars
1312  integer,parameter :: nqlwl0=0
1313  integer :: iq_ibz,isym_q,itim_q,npw,nomega,nqibz,mdf_type
1314  real(dp) :: eps_inf
1315  logical :: q_isirred
1316  !character(len=500) :: msg
1317 !arrays
1318  real(dp) :: qbz(3)
1319 
1320 ! *********************************************************************
1321 
1322  DBG_ENTER("COLL")
1323 
1324  npw = W%npw; nqibz = W%nqibz; nomega = W%nomega
1325 
1326  call get_bz_item(Qmesh,iq_bz,qbz,iq_ibz,isym_q,itim_q,isirred=q_isirred)
1327  !
1328  ! ========================================================
1329  ! ==== Branching for in-core or out-of-core solutions ====
1330  ! ========================================================
1331  if (screen_ihave_fgg(W,iq_ibz,how="Stored")) then
1332 
1333    if (q_isirred) then ! Symmetrization is not needed. Point the data in memory.
1334      call screen_fgg_qbz_set(W,iq_bz,nqlwl0,"Pointer")
1335      W%Fgg_qbz => W%Fgg(iq_ibz)
1336    else
1337      ! Allocate space. ! TODO Wings are not symmetrized but oh well
1338      call screen_fgg_qbz_set(W,iq_bz,nqlwl0,"Allocate")  ! Dimensions should not be changed.
1339 
1340      ! Out-of-place symmetrization.
1341      !em1_qibz => W%Fgg(iq_ibz)%mat; em1_qbz  => W%Fgg_qbz%mat
1342      call em1_symmetrize_op(iq_bz,npw,nomega,Gsph,Qmesh,W%Fgg(iq_ibz)%mat,W%Fgg_qbz%mat)
1343    end if
1344 
1345    if (W%has_ppmodel>0) then
1346      ! Symmetrize the ppmodel tables: em1_qibz => W%Fgg(iq_ibz)%mat
1347      call ppm_symmetrizer(W%PPm,iq_bz,Cryst,Qmesh,Gsph,npw,nomega,W%omega,W%Fgg(iq_ibz)%mat,&
1348 &      W%nfftf_tot,W%ngfftf,W%ae_rhor(:,1))
1349    end if
1350 
1351  else if (screen_ihave_fgg(W,iq_ibz,how="Allocated")) then
1352    MSG_ERROR("Fgg_iqibz is allocated but not initialized!")
1353 
1354  else
1355    if (W%fgg_qbz_idx /= iq_bz) then
1356      ! Must compute em1_qbz here. em1_qbz => W%Fgg_qbz%mat
1357 
1358      ! Allocate the BZ buffer.
1359      call screen_fgg_qbz_set(W,iq_bz,nqlwl0,"Allocate")
1360 
1361      if (W%Info%use_mdf /= MDL_NONE) then
1362        ! Compute the model-dielectric function at qbz on-the fly and in sequential
1363        !call wrtout(std_out,"Will compute MDF on the fly","COLL")
1364        call screen_mdielf(iq_bz,npw,nomega,W%Info%use_mdf,W%Info%eps_inf,Cryst,Qmesh,Vcp,Gsph,&
1365 &         W%nspden,W%nfftf_tot,W%ngfftf,W%ae_rhor,"EM1",W%Fgg_qbz%mat,xmpi_comm_self)
1366 
1367      else
1368        ! Read W(q_ibz) and symmetrize it (do this only if we don't have the correct q_bz in memory).
1369        call wrtout(std_out,sjoin("Out of core with file: ",W%fname),"COLL")
1370        call read_screening(em1_ncname,W%fname,npw,1,nomega,W%Fgg_qbz%mat,W%iomode,xmpi_comm_self,iqiA=iq_ibz)
1371 
1372        ! In-place symmetrization to get the q-point in the BZ.
1373        if (.not.q_isirred) then
1374          call em1_symmetrize_ip(iq_bz,npw,nomega,Gsph,Qmesh,W%Fgg_qbz%mat)
1375        end if
1376      end if
1377 
1378      W%Fgg_qbz%has_mat = MAT_STORED
1379    end if
1380 
1381    ABI_CHECK(W%Fgg_qbz%has_mat == MAT_STORED, "Wrong has_mat")
1382    !
1383    ! Ppmodel calculations with ppm tables in memory.
1384    ! TODO treat the case in which IBZ tables are stored in memory.
1385    if (W%has_ppmodel>0) then
1386      MSG_ERROR("Not implemented error")
1387      ! Symmetrize the ppmodel using em1_qibz.
1388      !%call ppm_symmetrizer(W%PPm,iq_bz,Cryst,Qmesh,Gsph,npw,nomega,W%omega,em1_qibz,W%nfftf_tot,W%ngfftf,W%ae_rhor(:,1))
1389      call ppm_symmetrizer(W%PPm,iq_bz,Cryst,Qmesh,Gsph,npw,nomega,W%omega,W%Fgg_qbz%mat ,W%nfftf_tot,W%ngfftf,W%ae_rhor(:,1))
1390    end if
1391  end if
1392 
1393  ! Calculate model dielectric function for this q-point in the BZ.
1394  eps_inf = W%Info%eps_inf; mdf_type = W%Info%use_mdf
1395 
1396  ! Model dielectric function. Only epsm-1 is supported here.
1397  !call wrtout(std_out," Calculating model dielectric function... ","COLL")
1398  !ABI_CHECK(W%nomega==1,"Cannot use nomega>1 in model dielectric function")
1399 
1400  ! W%Fgg_qbz%mat
1401 
1402  !%  call screen_mdielf(iq_bz,npw,nomega,mdf_type,eps_inf,Cryst,Qmesh,Vcp,Gsph,&
1403  !% &  W%nspden,W%nfftf_tot,W%ngfftf,W%ae_rhor,"EM1",em1_qbz,xmpi_comm_self)
1404 
1405  ! Store the index of the q-point in the BZ for checking purpose.
1406  W%fgg_qbz_idx = iq_bz
1407 
1408  DBG_EXIT("COLL")
1409 
1410 end subroutine screen_symmetrizer

m_screen/screen_t [ Types ]

[ Top ] [ m_screen ] [ Types ]

NAME

 screen_t

FUNCTION

  Object used to store the screening matrix in reciprocal space.

SOURCE

226  type,public :: screen_t
227 
228 ! scalars
229   integer :: iomode               ! Flag defining the IO mode.
230   integer :: debug_level=0        ! Internal Flag defining the debug level.
231   integer :: mqmem                ! =0 for out-of-core solution, =nqibz if entire matrix is stored in memory.
232   integer :: nI,nJ                ! Number of components (rows,columns) in chi|eps^-1. (1,1) if collinear.
233   integer :: nqibz                ! Number of q-points in the IBZ used.
234   integer :: nqlwl                ! Number of points used for the treatment of the long wave-length limit.
235   integer :: nomega               ! Total Number of frequencies used.
236   integer :: nomega_i             ! Number of purely imaginary frequencies used.
237   integer :: nomega_r             ! Number of real frequencies used.
238   integer :: npw                  ! Number of G vectors.
239   integer :: prtvol               ! Verbosity level.
240   integer :: has_ppmodel          ! 1 if PPmodel tables are stored.
241   integer :: has_fgg              ! 1 if Fgg tables are stored.
242   !%integer :: wmesh_type
243   !%real(dp) :: ecut
244   integer :: nfftf_tot
245   integer :: nspden
246 
247 !arrays
248   integer :: ngfftf(18)          ! Info on the FFT mesh used for ae_rhor (used for the model dielectric function)
249 
250   real(dp),allocatable :: ae_rhor(:,:)
251   ! ae_rhor(nfft,nspden)
252   ! Density in real space used to construct the TDDFT kernel or the model dielectric function.
253   ! NOTE that ae_rhor is given on the dense mesh as it contains the PAW onsite contribution.
254 
255   character(len=fnlen) :: fname=ABI_NOFILE  ! Name of the file used for the out-of-core solution.
256 
257   real(dp),allocatable :: qibz(:,:)
258   ! qibz(3,nqibz)
259   ! q-points in reduced coordinates
260 
261   !type(kmesh_t) :: Qmesh
262   ! Info the q-point sampling.
263 
264   real(dp),allocatable :: qlwl(:,:)
265   ! qlwl(3,nqlwl)
266   ! q-points used for the long wave-length limit treatment.
267 
268   complex(dpc),allocatable :: omega(:)
269   ! omega(nomega)
270   ! List of frequencies. Real frequencies are packed first.
271 
272   integer,allocatable :: gvec(:,:)
273   ! gvec(3,npw)
274   ! G-vectors used to describe the two-point function (r.l.u.).
275 
276   !type(gsphere_t) :: Gsphere
277   ! Info on the G-sphere. Note that the basis set does not depend on q.
278   ! The G-vectors are ordered in shells to facilitate the application of symmetry operations.
279   ! See m_gsphere.F90.
280 
281   logical,allocatable :: keep_q(:)
282    ! keep_q(nqibz)
283    ! Storage strategy: keep or not keep Em1(q) in memory.
284 
285   type(fgg_t),pointer :: Fgg(:)
286   ! Fgg(nqibz)
287   ! F_{G,G'}(q,w) for q in the IBZ.
288 
289   integer :: fgg_qbz_stat=FGG_QBZ_ISPOINTER
290   ! Status of Fgg_qbz
291 
292   integer :: fgg_qbz_idx=0
293   ! The index of the q-point in BZ pointed by Fgg_qbz. Used for debugging purpose.
294 
295   type(fgg_t),pointer :: Fgg_qbz  => null()
296   ! Buffer used for storing F_GG' at the point q_bz in the BZ
297   ! If q_bz is in the IBZ, Fgg_qbz *points* to Fgg(iq_ibz)
298   ! If q_bz is not in the IBZ, Fgg_qbz is *allocated* and used to store the symmetrized matrix.
299 
300   type(ppmodel_t) :: PPm
301   ! Structure storing the plasmon-pole parameters.
302 
303   type(screen_info_t) :: Info
304   ! Parameters used to construct the screening.
305 
306  end type screen_t
307 
308  public :: screen_nullify       ! Nullify all pointers before use.
309  public :: screen_init          ! Creation method.
310  !public :: screen_print         ! Print info
311  public :: screen_free          ! Free all associated pointers.
312  public :: screen_symmetrizer   ! Prepare the object for applying W_qbz.
313  public :: screen_w0gemv        ! Matrix vector multiplication \sum_{G'} F_{G,G') |u(G')>.
314  !public :: screen_times_ket     ! Compute \Sigma_c(\omega)|\phi> in reciprocal space.

m_screen/screen_w0gemv [ Functions ]

[ Top ] [ m_screen ] [ Functions ]

NAME

 screen_w0gemv

FUNCTION

  Perform the matrix multiplication W x vector in reciprocal space.

INPUTS

  W<screen_t>=
  in_npw=Number of G vectors in in_ket
  nspinor=Number of spinorial components.
  in_ket(in_npw)= |\phi> in reciprocal space.
  trans= On entry, TRANS specifies the operation to be performed as follows:
  TRANS = 'N' or 'n'   y := alpha*A*x + beta*y.
  TRANS = 'T' or 't'   y := alpha*A**T*x + beta*y.
  TRANS = 'C' or 'c'   y := alpha*A**H*x + beta*y.

OUTPUT

   out_ket(in_npw)= W |\phi\> in reciprocal space.
   ZGEMV  performs one of the matrix-vector operations
   *
   *     y := alpha*A*x + beta*y,   or   y := alpha*A**T*x + beta*y,   or
   *
   *     y := alpha*A**H*x + beta*y,
   *
   *  where alpha and beta are scalars, x and y are vectors and A is an m by n matrix.

PARENTS

      exc_build_block

CHILDREN

      get_bz_item,sqmat_itranspose

SOURCE

1450 subroutine screen_w0gemv(W,trans,in_npw,nspinor,only_diago,alpha,beta,in_ket,out_ket)
1451 
1452 
1453 !This section has been created automatically by the script Abilint (TD).
1454 !Do not modify the following lines by hand.
1455 #undef ABI_FUNC
1456 #define ABI_FUNC 'screen_w0gemv'
1457 !End of the abilint section
1458 
1459  implicit none
1460 
1461 !Arguments ------------------------------------
1462 !scalars
1463  integer,intent(in) :: in_npw,nspinor
1464  complex(gwpc),intent(in) :: alpha,beta
1465  logical,intent(in) :: only_diago
1466  character(len=*),intent(in) ::  trans
1467  type(screen_t),intent(in) :: W
1468 !arrays
1469  complex(gwpc),intent(in) :: in_ket(in_npw*nspinor)
1470  complex(gwpc),intent(out) :: out_ket(in_npw*nspinor)
1471 
1472 !Local variables-------------------------------
1473 !scalars
1474  integer :: ig,lda
1475 !arrays
1476  complex(gwpc),ABI_CONTIGUOUS pointer :: em1_qbz(:,:)
1477 
1478 ! *************************************************************************
1479 
1480  lda = W%npw
1481  em1_qbz => W%Fgg_qbz%mat(:,:,1)
1482 
1483  if (.not.only_diago) then
1484    !out_ket = MATMUL(TRANSPOSE(CONJG(em1_qbz(1:in_npw,1:in_npw)),in_ket)
1485    call xgemv(trans,in_npw,in_npw,alpha,em1_qbz,lda,in_ket,1,beta,out_ket,1)
1486 
1487  else
1488    if (beta /= czero_gw) then
1489      !
1490      if ( firstchar(trans,(/"C"/)) ) then
1491        do ig=1,in_npw
1492          out_ket(ig) = alpha * CONJG(em1_qbz(ig,ig)) * in_ket(ig) + beta * out_ket(ig)
1493        end do
1494      else if ( firstchar(trans,(/"N","T"/)) ) then
1495        do ig=1,in_npw
1496          out_ket(ig) = alpha * em1_qbz(ig,ig) * in_ket(ig) + beta * out_ket(ig)
1497        end do
1498      else
1499        MSG_ERROR(sjoin("Wrong trans:", trans))
1500      end if
1501 
1502    else
1503      ! beta==0
1504      if ( firstchar(trans,(/"C"/)) ) then
1505        do ig=1,in_npw
1506          out_ket(ig) = alpha * CONJG(em1_qbz(ig,ig)) * in_ket(ig)
1507        end do
1508      else if ( firstchar(trans,(/"N","T"/)) ) then
1509        do ig=1,in_npw
1510          out_ket(ig) = alpha * em1_qbz(ig,ig) * in_ket(ig)
1511        end do
1512      else
1513        MSG_ERROR(sjoin("Wrong trans:", trans))
1514      end if
1515 
1516    end if
1517  end if
1518 
1519 end subroutine screen_w0gemv