TABLE OF CONTENTS


ABINIT/lin_interpq_gam [ Functions ]

[ Top ] [ Functions ]

NAME

 lin_interpq_gam

FUNCTION

 This routine interpolates the electron phonon coupling matrix
 to a give _qpoint_ by linear methods

COPYRIGHT

 Copyright (C) 2010-2018 ABINIT group (MVer)
 This file is distributed under the terms of the
 GNU General Public Licence, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

INPUTS

   nbranch=Number of phonon branches 3*natom
   nqbz=Number of q-points in full BZ.
   gamma_qpt(2,nbranch**2,nsppol,nqbz)=Gamma matrix in full BZ.
   nsppol=Number of independent spin polarizations.
   qpt = q-point in reciprocal space for interpolation

OUTPUT

   gam_now = interpolated gamma matrix at qpt

PARENTS

CHILDREN

      wrap2_zero_one

SOURCE

35 #if defined HAVE_CONFIG_H
36 #include "config.h"
37 #endif
38 
39 #include "abi_common.h"
40 
41 subroutine lin_interpq_gam(gamma_qpt,nbranch,nqbz,nsppol,gam_now,isppol,kptrlatt,qpt)
42 
43  use defs_basis
44  use defs_elphon
45  use m_errors
46  use m_profiling_abi
47 
48  use m_numeric_tools,  only : wrap2_zero_one, interpol3d
49 
50 !This section has been created automatically by the script Abilint (TD).
51 !Do not modify the following lines by hand.
52 #undef ABI_FUNC
53 #define ABI_FUNC 'lin_interpq_gam'
54 !End of the abilint section
55 
56  implicit none
57 
58 !Arguments ------------------------------------
59 !scalars
60  integer, intent(in) :: nbranch,isppol,nsppol,nqbz
61 !arrays
62  integer, intent(in) :: kptrlatt(3,3)
63  real(dp), intent(in) :: qpt(3)
64  real(dp),intent(in) :: gamma_qpt(2,nbranch**2,nsppol,nqbz)
65  real(dp), intent(out) :: gam_now(2,nbranch**2)
66 
67 !Local variables-------------------------------
68 !scalars
69  integer :: ibranch
70  real(dp) :: res
71  character(len=500) :: msg
72 !arrays
73  real(dp) :: qpt_incube(3)
74 
75 ! *************************************************************************
76 
77 !the array gamma_qpt has dimensions: (2,nbranch**2,nsppol,elph_ds%k_fine%nkpt)
78  if (kptrlatt(1,1)*kptrlatt(2,2)*kptrlatt(3,3) /= nqbz) then
79    write(msg,'(a,2i0)')"Wrong dimensions in gamma_qpt ",kptrlatt(1,1)*kptrlatt(2,2)*kptrlatt(3,3),nqbz
80    MSG_ERROR(msg)
81  end if
82 
83  call wrap2_zero_one(qpt(1),qpt_incube(1),res)
84  call wrap2_zero_one(qpt(2),qpt_incube(2),res)
85  call wrap2_zero_one(qpt(3),qpt_incube(3),res)
86  
87  do ibranch=1,nbranch**2
88    gam_now(1,ibranch) = interpol3d(qpt_incube,kptrlatt(1,1),kptrlatt(2,2),kptrlatt(3,3),gamma_qpt(1,ibranch,isppol,:))
89    gam_now(2,ibranch) = interpol3d(qpt_incube,kptrlatt(1,1),kptrlatt(2,2),kptrlatt(3,3),gamma_qpt(2,ibranch,isppol,:))
90  end do
91 
92 end subroutine lin_interpq_gam