Commit 7e048161 authored by Samuel Poncé's avatar Samuel Poncé

New EPW feature: calculation of carrier mobilities.

parent 6ad84919
......@@ -19,7 +19,7 @@ MODFLAGS = -I../../iotk/src -I../../UtilXlib/ -I../../Modules -I../../KS_Solvers
EPWOBJS = \
transportcom.o \
io_epw.o\
io_epw.o \
elph2.o \
a2f.o \
allocate_epwq.o \
......@@ -77,6 +77,7 @@ kernels_iso_iaxis.o \
kernels_raxis.o \
kfold.o \
kpointdivision.o \
kpoint_grid_epw.o \
ktokpmq.o \
loadkmesh.o \
loadqmesh.o \
......@@ -112,12 +113,15 @@ star_q2.o \
stop_epw.o \
vmebloch2wan.o \
vmewan2bloch.o \
wannier.o \
wannierEPW.o \
wannierize.o \
pw2wan90epw.o \
wigner_seitz2.o \
wigner_seitz.o \
write_ephmat.o \
iterativebte.o \
transport_coeffs.o \
scattering_rate.o \
io_scattering.o \
system_mem_usage.o
......
......@@ -52,6 +52,7 @@
!
! ... zero up to a given accuracy
!
REAL(DP), PARAMETER :: eps2 = 1.0E-2_DP
REAL(DP), PARAMETER :: eps4 = 1.0E-4_DP
REAL(DP), PARAMETER :: eps6 = 1.0E-6_DP
REAL(DP), PARAMETER :: eps8 = 1.0E-8_DP
......
......@@ -38,6 +38,7 @@
USE modes, ONLY : t, npert, u, name_rap_mode, num_rap_mode
USE qpoint, ONLY : eigqts, igkq
USE klist, ONLY : nks
USE transportcom, ONLY : transp_temp
!
IMPLICIT NONE
INTEGER :: ik, ipol
......@@ -59,6 +60,7 @@
IF(ALLOCATED(xk_all)) DEALLOCATE (xk_all)
IF(ALLOCATED(et_all)) DEALLOCATE (et_all)
IF(ALLOCATED(wslen)) DEALLOCATE (wslen)
IF(ALLOCATED(eps_rpa)) DEALLOCATE (eps_rpa)
IF(ALLOCATED(eps_rpa)) DEALLOCATE (eps_rpa)
!
ELSE
......
......@@ -60,6 +60,9 @@
sigmar_all(:,:), &! To store sigmar, sigmai and zi globally
sigmai_all(:,:), &!
sigmai_mode(:,:,:), &!
Fi_all(:,:,:), &!
F_current(:,:,:), &!
F_SERTA(:,:,:), &!
zi_all(:,:), &!
esigmar_all(:,:,:), &!
esigmai_all(:,:,:), &!
......@@ -71,6 +74,9 @@
zstar(:,:,:), &! Born effective charges
epsi(:,:), &! dielectric tensor
inv_tau_all(:,:,:), &! scattering rate
inv_tau_allcb(:,:,:), &! Second scattering rate (for both)
zi_allvb(:,:,:), &! Z-factor in scattering rate
zi_allcb(:,:,:), &! Second Z-factor in scattering rate (for both VB and CB calculations)
ifc(:,:,:,:,:,:,:) ! Interatomic force constant in real space
REAL(KIND=DP) :: &!
efnew ! SP: Fermi level on the fine grid. Added globaly for efficiency reason
......
......@@ -105,9 +105,9 @@
!! true if we are using time reversal
!
! Local variables
logical :: exst!, addnlcc
LOGICAL :: exst!, addnlcc
!
integer :: ik, ipert, mode, ibnd, jbnd, ig, nkq, ipool, &
INTEGER :: ik, ipert, mode, ibnd, jbnd, ig, nkq, ipool, &
ik0, igkq_tmp (npwx), imap, &
ipooltmp, nkq_abs, ipol, npw
INTEGER :: ng0vec, ngxx, lower_bnd, upper_bnd, nkk, nkk_abs
......@@ -115,7 +115,7 @@
! bound for the allocation of the array gmap
!
! variables for rotating the wavefunctions (in order to USE q in the irr wedge)
real(kind=DP) :: xktmp(3), sxk(3)
REAL(kind=DP) :: xktmp(3), sxk(3)
REAL(kind=DP) :: g0vec_all_r(3,125)
!! Variables for folding of k+q grid
REAL(kind=DP) :: zero_vect(3)
......
This diff is collapsed.
......@@ -55,9 +55,9 @@
pwc, nswc, nswfc, nswi, filukq, filukk, &
nbndsub, nbndskip, system_2d, delta_approx, &
title, int_mob, scissor, iterative_bte, scattering, &
ncarrier, carrier, scattering_serta, &
scattering_0rta, longrange, shortrange,restart, &
restart_freq, prtgkk, nel, meff, epsiHEG
ncarrier, carrier, scattering_serta, restart, restart_freq, &
scattering_0rta, longrange, shortrange, scatread, &
restart_filq, prtgkk, nel, meff, epsiHEG
USE elph2, ONLY : elph
USE start_k, ONLY : nk1, nk2, nk3
USE constants_epw, ONLY : ryd2mev, ryd2ev, ev2cmm1, kelvin2eV
......@@ -121,7 +121,7 @@
specfun_el, specfun_ph, wmin_specfun, wmax_specfun, nw_specfun, &
delta_approx, scattering, int_mob, scissor, ncarrier, carrier, &
iterative_bte, scattering_serta, scattering_0rta, longrange, shortrange,&
restart, restart_freq, prtgkk, nel, meff, epsiHEG
scatread, restart, restart_freq, restart_filq, prtgkk, nel, meff, epsiHEG
! tphases, fildvscf0
!
......@@ -252,8 +252,10 @@
! specfun_pl : if .TRUE. calculate plason spectral function
! restart : if .true. a run can be restarted from the interpolation level
! restart_freq : Create a restart point every restart_freq q/k-points
! restart_filq : Use to merge different q-grid scattering rates (name of the file)
! scattering : if .true. scattering rates are calculated
! scattering_serta: if .true. scattering rates are calculated using self-energy relaxation-time-approx
! scatread : if .true. the current scattering rate file is read from file.
! scattering_0rta : if .true. scattering rates are calculated using 0th order relaxation-time-approx
! int_mob : if .true. computes the intrinsic mobilities. This means that the
! electron and hole carrier density is equal.
......@@ -362,6 +364,7 @@
dvscf_dir = '.'
prefix = 'pwscf'
filqf = ' '
restart_filq = ' '
filkf = ' '
fildrho = ' '
fildvscf = ' '
......@@ -449,6 +452,7 @@
system_2d = .false.
scattering = .false.
scattering_serta = .false.
scatread = .false.
scattering_0rta = .false.
int_mob = .false.
iterative_bte = .false.
......@@ -570,12 +574,14 @@
CALL errore('epw_init', 'tempsmax should be greater than tempsmin',1)
IF ( int_mob .AND. efermi_read) CALL errore('epw_init', &
'Fermi level can not be set (efermi_read) when computing intrinsic mobilities',1)
IF ( int_mob .AND. (ABS(ncarrier) > 1E+5) ) CALL errore('epw_init', &
'You cannot compute intrinsic mobilities and doped mobilities at the same time',1)
! IF ( int_mob .AND. (ABS(ncarrier) > 1E+5) ) CALL errore('epw_init', &
! 'You cannot compute intrinsic mobilities and doped mobilities at the same time',1)
IF ( (ABS(ncarrier) > 1E+5) .and. .not. carrier ) CALL errore('epw_init', &
'carrier must be .true. if you specify ncarrier.',1)
IF ( carrier .AND. (ABS(ncarrier) < 1E+5) ) CALL errore('epw_init', &
'The absolute value of the doping carrier concentration must be larger than 1E5 cm^-3',1)
IF ( (iterative_bte .AND. scattering_serta) .OR. (iterative_bte .AND.scattering_0rta) ) CALL errore('epw_init', &
'You should first do a run in the SRTA to get the initial scattering_rate files.',1)
! IF ( (iterative_bte .AND. scattering_serta) .OR. (iterative_bte .AND.scattering_0rta) ) CALL errore('epw_init', &
! 'You should first do a run in the SRTA to get the initial scattering_rate files.',1)
IF ( (longrange .OR. shortrange) .AND. (.not. lpolar)) CALL errore('epw_init',&
&'Error: longrange or shortrange can only be true if lpolar is true as well.',1)
IF ( longrange .AND. shortrange) CALL errore('epw_init',&
......
......@@ -30,6 +30,7 @@
USE uspp_param, ONLY : upf
USE spin_orb, ONLY : domag
USE constants, ONLY : degspin, pi
USE constants_epw, ONLY : zero
USE noncollin_module, ONLY : noncolin, m_loc, angle1, angle2, ux
USE wvfct, ONLY : nbnd, et
USE output, ONLY : fildrho
......@@ -42,8 +43,7 @@
USE control_lr, ONLY : alpha_pv, nbnd_occ
USE gamma_gamma, ONLY : has_equivalent, asr, nasr, n_diff_sites, &
equiv_atoms, n_equiv_atoms, with_symmetry
USE partial, ONLY : &
done_irr
USE partial, ONLY : done_irr
USE modes, ONLY : u, npertx, npert, nirr, nmodes, num_rap_mode
USE lr_symm_base, ONLY : gi, gimq, irotmq, minus_q, nsymq, invsymq, rtau
USE qpoint, ONLY : xq
......@@ -53,10 +53,12 @@
USE mp, ONLY : mp_bcast
USE mp, ONLY : mp_max, mp_min
USE mp_pools, ONLY : inter_pool_comm
USE epwcom, ONLY : xk_cryst
USE epwcom, ONLY : xk_cryst,scattering, nstemp, tempsmin, tempsmax, &
temps
USE fft_base, ONLY : dfftp
USE gvecs, ONLY : doublegrid
USE start_k, ONLY : nk1, nk2, nk3
USE transportcom, ONLY : transp_temp
!
implicit none
!
......@@ -91,6 +93,8 @@
INTEGER :: js
!! counter on atomic type
INTEGER :: last_irr_eff
!! Last effective irr
INTEGER :: itemp
!!
REAL(kind=DP) :: rhotot
!! total charge
......@@ -392,7 +396,27 @@
DO irr = 1, nirr
npertx = max (npertx, npert (irr) )
ENDDO
!
transp_temp(:) = zero
! In case of scattering calculation
IF ( scattering ) THEN
!
IF ( maxval(temps(:)) > zero ) THEN
transp_temp(:) = temps(:)
ELSE
IF ( nstemp .eq. 1 ) THEN
transp_temp(1) = tempsmin
ELSE
DO itemp = 1, nstemp
transp_temp(itemp) = tempsmin + dble(itemp-1) * &
( tempsmax - tempsmin ) / dble(nstemp-1)
ENDDO
ENDIF
ENDIF
ENDIF
! We have to bcast here because before it has not been allocated
CALL mp_bcast (transp_temp, ionode_id, world_comm) !
!
CALL stop_clock ('epw_setup')
RETURN
!
......@@ -448,6 +472,3 @@
!-----------------------------------------------------------------------
END SUBROUTINE epw_setup_restart
!-----------------------------------------------------------------------
......@@ -266,6 +266,8 @@
!! if .true. scattering rates are calculated
LOGICAL :: scattering_serta
!! if .true. scattering rates are calculated using self-energy relaxation-time-approx
LOGICAL :: scatread
!! if .true. the scattering rates are read from file.
LOGICAL :: scattering_0rta
!! if .true. scattering rates are calculated using 0th order relaxation-time-approx
LOGICAL :: int_mob
......@@ -343,6 +345,8 @@ MODULE output_epw
!! output file for the deltavscf used as a fake perturbation to set phases
CHARACTER (LEN=80) :: fila2f
!! input file containing eliashberg spectral function
CHARACTER (LEN=80) :: restart_filq
!! input file to restart from an exisiting q-file
!
END MODULE output_epw
!
......
This diff is collapsed.
This diff is collapsed.
!
! Copyright (C) 2010-2016 Samuel Ponce', Roxana Margine, Carla Verdi, Feliciano Giustino
!
! This file is distributed under the terms of the GNU General Public
! License. See the file `LICENSE' in the root directory of the
! present distribution, or http://www.gnu.org/copyleft.gpl.txt .
!
!-----------------------------------------------------------------------
SUBROUTINE kpoint_grid_epw ( nrot, time_reversal, skip_equivalence, s, t_rev, &
bg, npk, nk1, nk2, nk3, nks, xk, wk, BZtoIBZ, s_BZtoIBZ)
!-----------------------------------------------------------------------
!!
!! Automatic generation of a uniform grid of k-points with symmetry.
!! Routine copied from PW/src/kpoint_grid.f90.
!! We had to duplicate because the BZtoIBZ array was deallocated and is needed in
!! EPW
!!
USE kinds, ONLY: DP
IMPLICIT NONE
!
INTEGER, INTENT(in) :: nrot, npk, nk1, nk2, nk3, &
t_rev(48), s(3,3,48)
INTEGER, INTENT(inout) :: s_BZtoIBZ(3,3,nk1*nk2*nk3)
!! Symeetry matrix that links an point to its IBZ friend.
INTEGER, INTENT(inout) :: BZtoIBZ (nk1*nk2*nk3)
!! Number of rotation
LOGICAL, INTENT(in) :: time_reversal
!! True if time reversal
LOGICAL, INTENT(in) :: skip_equivalence
!! True if equivalent point
REAL(DP), INTENT(in) :: bg(3,3)
!! Reciprocal space vectors
!
INTEGER, INTENT(out) :: nks
real(DP), INTENT(out):: xk(3,npk)
real(DP), INTENT(out):: wk(npk)
! LOCAL:
INTEGER :: s_save(3,3,nk1*nk2*nk3)
real(DP), PARAMETER :: eps=1.0d-5
real(DP) :: xkr(3), fact, xx, yy, zz
real(DP), ALLOCATABLE:: xkg(:,:), wkk(:)
INTEGER :: nkr, i,j,k, ns, n, nk, equiv(nk1*nk2*nk3), ik
LOGICAL :: in_the_list
!
nkr=nk1*nk2*nk3
ALLOCATE (xkg( 3,nkr),wkk(nkr))
equiv(:) = 0
s_save(:,:,:) = 0
!
DO i=1,nk1
DO j=1,nk2
DO k=1,nk3
! this is nothing but consecutive ordering
n = (k-1) + (j-1)*nk3 + (i-1)*nk2*nk3 + 1
! xkg are the components of the complete grid in crystal axis
xkg(1,n) = dble(i-1)/nk1
xkg(2,n) = dble(j-1)/nk2
xkg(3,n) = dble(k-1)/nk3
ENDDO
ENDDO
ENDDO
! equiv(nk) =nk : k-point nk is not equivalent to any previous k-point
! equiv(nk)!=nk : k-point nk is equivalent to k-point equiv(nk)
DO nk=1,nkr
equiv(nk)=nk
ENDDO
IF ( skip_equivalence ) THEN
CALL infomsg('kpoint_grid', 'ATTENTION: skip check of k-points equivalence')
wkk = 1.d0
ELSE
DO nk=1,nkr
! check if this k-point has already been found equivalent to another
IF (equiv(nk) == nk) THEN
wkk(nk) = 1.0d0
! check if there are equivalent k-point to this in the list
! (excepted those previously found to be equivalent to another)
! check both k and -k
DO ns=1,nrot
DO i=1,3
xkr(i) = s(i,1,ns) * xkg(1,nk) &
+ s(i,2,ns) * xkg(2,nk) &
+ s(i,3,ns) * xkg(3,nk)
xkr(i) = xkr(i) - nint( xkr(i) )
ENDDO
IF(t_rev(ns)==1) xkr = -xkr
xx = xkr(1)*nk1
yy = xkr(2)*nk2
zz = xkr(3)*nk3
in_the_list = abs(xx-nint(xx))<=eps .and. &
abs(yy-nint(yy))<=eps .and. &
abs(zz-nint(zz))<=eps
IF (in_the_list) THEN
i = mod ( nint ( xkr(1)*nk1 + 2*nk1), nk1 ) + 1
j = mod ( nint ( xkr(2)*nk2 + 2*nk2), nk2 ) + 1
k = mod ( nint ( xkr(3)*nk3 + 2*nk3), nk3 ) + 1
n = (k-1) + (j-1)*nk3 + (i-1)*nk2*nk3 + 1
IF (n>nk .and. equiv(n)==n) THEN
equiv(n) = nk
wkk(nk)=wkk(nk)+1.0d0
s_save(:,:,n) = s(:,:,ns)
ELSE
IF (equiv(n)/=nk .or. n<nk ) CALL errore('kpoint_grid', &
'something wrong in the checking algorithm',1)
ENDIF
ENDIF
IF ( time_reversal ) THEN
xx =-xkr(1)*nk1
yy =-xkr(2)*nk2
zz =-xkr(3)*nk3
in_the_list=abs(xx-nint(xx))<=eps.and.abs(yy-nint(yy))<=eps &
.and. abs(zz-nint(zz))<=eps
IF (in_the_list) THEN
i = mod ( nint (-xkr(1)*nk1 + 2*nk1), nk1 ) + 1
j = mod ( nint (-xkr(2)*nk2 + 2*nk2), nk2 ) + 1
k = mod ( nint (-xkr(3)*nk3 + 2*nk3), nk3 ) + 1
n = (k-1) + (j-1)*nk3 + (i-1)*nk2*nk3 + 1
IF (n>nk .and. equiv(n)==n) THEN
equiv(n) = nk
wkk(nk)=wkk(nk)+1.0d0
s_save(:,:,n) = -s(:,:,ns)
ELSE
IF (equiv(n)/=nk.or.n<nk) CALL errore('kpoint_grid', &
'something wrong in the checking algorithm',2)
ENDIF
ENDIF
ENDIF
ENDDO
ENDIF
ENDDO
ENDIF
! count irreducible points and order them
nks=0
fact=0.0d0
!
DO nk=1,nkr
BZtoIBZ(nk) = equiv(nk)
ENDDO
!
DO nk=1,nkr
IF (equiv(nk) == nk) THEN
nks=nks+1
IF (nks>npk) CALL errore('kpoint_grid','too many k-points',1)
wk(nks) = wkk(nk)
fact = fact+wk(nks)
! bring back into to the first BZ
DO i=1,3
xk(i,nks) = xkg(i,nk)-nint(xkg(i,nk))
ENDDO
! DBSP
BZtoIBZ(nk) = nks
! Change all the one above
DO ik=nk,nkr
IF (equiv(ik) == nk) THEN
BZtoIBZ(ik) = nks
ENDIF
ENDDO
ENDIF
ENDDO
! Now do the symmetry mapping.
DO nk=1,nkr
! If its an irreducible point
IF ( equiv(nk) == nk ) THEN
! Then you have the identity matrix
s_BZtoIBZ(:,:,nk) = s(:,:,1)
ELSE
s_BZtoIBZ(:,:,nk) = s_save(:,:,nk)
ENDIF
ENDDO
!
! Store mapping from all the point to IBZ
! go to cartesian axis (in units 2pi/a0)
CALL cryst_to_cart(nks,xk,bg,1)
! normalize weights to one
DO nk=1,nks
wk(nk) = wk(nk)/fact
ENDDO
DEALLOCATE(xkg,wkk)
RETURN
END SUBROUTINE kpoint_grid_epw
......@@ -161,16 +161,208 @@
!
nkq_abs = n
!
end subroutine ktokpmq
!--------------------------------------------------------
END SUBROUTINE ktokpmq
!--------------------------------------------------------
!
!--------------------------------------------------------
SUBROUTINE ktokpmq_fine (xkf_cryst, xk, xq, sign, ipool, nkq, nkq_abs)
!--------------------------------------------------------
!!
!! For a given k point in cart coord, find the index
!! of the corresponding (k + sign*q) point on the fine
!! homogeneous grid
!!
!! In the parallel case, determine also the pool number
!! nkq is the in-pool index, nkq_abs is the absolute
!! index
!!
!--------------------------------------------------------
!
USE kinds, ONLY : DP
USE cell_base, ONLY : at
USE epwcom, ONLY : nkf1, nkf2, nkf3
USE elph2, ONLY : nkqtotf
USE mp_global, ONLY : nproc_pool, npool
USE mp_images, ONLY : nproc_image
USE mp, ONLY : mp_barrier, mp_bcast
!
USE mp_world, ONLY : mpime
! m
implicit none
!
INTEGER, INTENT (in) :: sign
!! +1 for searching k+q, -1 for k-q
INTEGER, INTENT (out) :: nkq
!! in the parallel case, the pool hosting the k+-q point
INTEGER, INTENT (out) :: nkq_abs
!! the index of k+sign*q
INTEGER, INTENT (out) :: ipool
!! Index of the pool where the point is
!
REAL(kind=DP), INTENT (in) :: xkf_cryst(3,nkqtotf/2)
!! coordinates of k points and q point
REAL(kind=DP), INTENT (in) :: xk(3)
!! coordinates of k points and q points
REAL(kind=DP), INTENT (in) :: xq(3)
!! Coordinates of k+q point
!
! work variables
!
INTEGER :: kunit
!! the absolute index of k+sign*q (in the full k grid)
real(kind=DP) :: xxk (3), xxq (3)
integer :: n, ik
real(kind=DP) :: xx, yy, zz, xx_c, yy_c, zz_c, eps
logical :: in_the_list, found
!
integer :: iks, nkl, nkr, jpool
!
kunit =1
! loosy tolerance, no problem since we use integer comparisons
eps = 1.d-5
IF (abs(sign).ne.1) call errore('ktokpmq_fine','sign must be +1 or -1',1)
!
! bring k and q in crystal coordinates
!
xxk = xk
xxq = xq
!
CALL cryst_to_cart (1, xxk, at, -1)
CALL cryst_to_cart (1, xxq, at, -1)
!
! check that k is actually on a uniform mesh centered at gamma
!
!print*,'xxk ',xxk(:)
!IF (mpime == 0) THEN
! write(910,*)'xxk ',xxk
! write(910,*)'xxq ',xxq
! write(910,*)'nkf1 ',nkf1
!ENDIf
!print*,'xxk ',xxk(:)
!IF (mpime == 1) THEN
! write(911,*)'xxk ',xxk
! write(911,*)'xxq ',xxq
! write(911,*)'nkf1 ',nkf1
! write(911,*)'xkf_cryst ',xkf_cryst(:,1)
! write(911,*)'xkf_cryst ',xkf_cryst(:,2)
! write(911,*)'xkf_cryst ',xkf_cryst(:,3)
! write(911,*)'xkf_cryst ',xkf_cryst(:,4)
!ENDIF
!---------------------------------
subroutine ckbounds(lower, upper)
!---------------------------------
!!
!! Subroutine finds the lower and upper
!! bounds of the coarse k-grid in parallel
!!
!---------------------------------
xx = xxk(1)*nkf1
yy = xxk(2)*nkf2
zz = xxk(3)*nkf3
in_the_list = abs(xx-nint(xx)).le.eps .and. &
abs(yy-nint(yy)).le.eps .and. &
abs(zz-nint(zz)).le.eps
IF (.not.in_the_list) call errore('ktokpmq_fine','is this a uniform k-mesh?',1)
!
IF ( xx .lt. -eps .or. yy .lt. -eps .or. zz .lt. -eps ) &
call errore('ktokpmq_fine','coarse k-mesh needs to be strictly positive in 1st BZ',1)
!
! now add the phonon wavevector and check that k+q falls again on the k grid
!
xxk = xxk + dble(sign) * xxq
!
xx = xxk(1)*nkf1
yy = xxk(2)*nkf2
zz = xxk(3)*nkf3
in_the_list = abs(xx-nint(xx)).le.eps .and. &
abs(yy-nint(yy)).le.eps .and. &
abs(zz-nint(zz)).le.eps
IF (.not.in_the_list) call errore('ktokpmq_fine','k+q does not fall on k-grid',1)
!
! find the index of this k+q in the k-grid
!
! make sure xx, yy and zz are in 1st BZ
!
CALL backtoBZ( xx, yy, zz, nkf1, nkf2, nkf3 )
!
n = 0
found = .false.
DO ik = 1, nkqtotf/2
xx_c = xkf_cryst(1,ik)*nkf1
yy_c = xkf_cryst(2,ik)*nkf2
zz_c = xkf_cryst(3,ik)*nkf3
!
! IF (mpime == 1) THEN
! write(911,*)'ik ',ik
! write(911,*)'xx_c ',xx_c
! write(911,*)'yy_c ',yy_c
! write(911,*)'zz_c ',zz_c
! ENDIF
!
! check that the k-mesh was defined in the positive region of 1st BZ
!
IF ( xx_c .lt. -eps .or. yy_c .lt. -eps .or. zz_c .lt. -eps ) &
call errore('ktokpmq_fine','coarse k-mesh needs to be strictly positive in 1st BZ',1)
!
found = nint(xx_c) .eq. nint(xx) .and. &
nint(yy_c) .eq. nint(yy) .and. &
nint(zz_c) .eq. nint(zz)
IF (found) THEN
n = ik
EXIT
ENDIF
ENDDO
!
! 26/06/2012 RM
! since coarse k- and q- meshes are commensurate, one can easily find n
!
! n = nint(xx) * nk2 * nk3 + nint(yy) * nk3 + nint(zz) + 1
!
IF (n .eq. 0) call errore('ktokpmq_fine','problem indexing k+q',1)
!
! Now n represents the index of k+sign*q in the original k grid.
! In the parallel case we have to find the corresponding pool
! and index in the pool
!
#ifdef __MPI
!
npool = nproc_image/nproc_pool
!
DO jpool = 0, npool-1
!
nkl = kunit * ( nkqtotf / npool )
nkr = ( nkqtotf - nkl * npool ) / kunit
!
! the reminder goes to the first nkr pools (0...nkr-1)
!
IF ( jpool < nkr ) nkl = nkl + kunit
!
! the index of the first k point in this pool
!
iks = nkl * jpool + 1
IF ( jpool >= nkr ) iks = iks + nkr * kunit
!
IF (n.ge.iks) then
ipool = jpool+1
nkq = n - iks + 1
ENDIF
!
ENDDO
!
#else
!
nkq = n
!
#endif
!
nkq_abs = n
!
!--------------------------------------------------------------
END SUBROUTINE ktokpmq_fine
!--------------------------------------------------------------
!
!---------------------------------
SUBROUTINE ckbounds(lower, upper)
!---------------------------------
!!
!! Subroutine finds the lower and upper
!! bounds of the coarse k-grid in parallel
!!
!---------------------------------
!
use pwcom, ONLY : nkstot
use mp_global, ONLY : my_pool_id, npool
......@@ -198,23 +390,30 @@ subroutine ckbounds(lower, upper)
lower = 1
upper = nkstot
#endif
end subroutine
!---------------------------------
subroutine para_bounds(lower, upper, total)
!---------------------------------
!!
!! Subroutine finds the lower and upper
!! bounds if we split some quantity over pools
!!
!---------------------------------
!
use mp_global, ONLY : my_pool_id, npool
!-------------------------------------
END SUBROUTINE ckbounds
!--------------------------------------
!
!---------------------------------