TABLE OF CONTENTS


ABINIT/m_prep_kgb [ Modules ]

[ Top ] [ Modules ]

NAME

  m_prep_kgb

FUNCTION

COPYRIGHT

  Copyright (C) 1998-2018 ABINIT group (FBottin,MT,GZ,MD,FDahm)
  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

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_prep_kgb
28 
29  use defs_basis
30  use defs_abitypes
31  use m_abicore
32  use m_errors
33  use m_xmpi
34 
35  use m_time,        only : timab
36  use m_bandfft_kpt, only : bandfft_kpt, bandfft_kpt_get_ikpt, bandfft_kpt_type
37  use m_pawcprj,     only : pawcprj_type
38  use m_hamiltonian, only : gs_hamiltonian_type
39  use m_nonlop,      only : nonlop
40  use m_getghc,      only : multithreaded_getghc
41  use m_fft,         only : fourwf
42 
43  implicit none
44 
45  private

ABINIT/prep_fourwf [ Functions ]

[ Top ] [ Functions ]

NAME

 prep_fourwf

FUNCTION

 this routine prepares the data to the call of fourwf.

INPUTS

  blocksize= size of block for FFT
  cwavef(2,npw*ndat)=planewave coefficients of wavefunction (one spinorial component?).
  dtfil <type(datafiles_type)>=variables related to files
  gvnlc=matrix elements <G|Vnonlocal|C>
  kg_k(3,npw_k)=reduced planewave coordinates.
  lmnmax=if useylm=1, max number of (l,m,n) comp. over all type of psps
        =if useylm=0, max number of (l,n)   comp. over all type of psps
  mgfft=maximum size of 1d ffts
  mpi_enreg=informations about mpi parallelization
  mpsang= 1+maximum angular momentum for nonlocal pseudopotentials
  mpssoang= 1+maximum (spin*angular momentum) for nonlocal pseudopotentials
  natom=number of atoms in cell.
  nband_k=number of bands at this k point for that spin polarization
  ndat=number of FFT to do in //
  ngfft(18)= contain all needed information about 3D FFT
  npw_k=number of plane waves at this k point
  nspinor=number of spinorial components of the wavefunctions
  ntypat=number of types of atoms in unit cell.
  n4,n5,n6 used for dimensionning of vlocal
  option_fourwf=option for fourwf (see fourwf.F90)
  prtvol=control print volume and debugging output
  ucvol=unit cell volume
  [bandfft_kpt_tab]= (optional) if present, contains tabs used to implement
                     the "band-fft" parallelism
                      if not present, the bandfft_kpt global variable is used
  [use_gpu_cuda]= (optional) 1 if Cuda (GPU) is on

OUTPUT

  gwavef=(2,npw*ndat)=matrix elements <G|H|C>.

SIDE EFFECTS

PARENTS

      mkrho,posdoppler,vtowfk

CHILDREN

      fourwf,gpu_fourwf,prep_index_wavef_bandpp,prep_wavef_sym_do
      prep_wavef_sym_undo,timab,xmpi_alltoallv,xmpi_sum

SOURCE

 929 subroutine prep_fourwf(rhoaug,blocksize,cwavef,wfraug,iblock,istwf_k,mgfft,&
 930 &          mpi_enreg,nband_k,ndat,ngfft,npw_k,n4,n5,n6,occ_k,option_fourwf,ucvol,wtk,&
 931 &          bandfft_kpt_tab,use_gpu_cuda) ! Optional arguments
 932 
 933 
 934 !This section has been created automatically by the script Abilint (TD).
 935 !Do not modify the following lines by hand.
 936 #undef ABI_FUNC
 937 #define ABI_FUNC 'prep_fourwf'
 938 !End of the abilint section
 939 
 940  implicit none
 941 
 942 !Arguments ------------------------------------
 943 !scalars
 944  integer,intent(in) :: blocksize,iblock,istwf_k,mgfft,n4,n5,n6,nband_k,ndat,npw_k
 945  integer,intent(in) :: option_fourwf
 946  integer,intent(in),optional :: use_gpu_cuda
 947  real(dp),intent(in) :: ucvol,wtk
 948  type(bandfft_kpt_type),optional,target,intent(in) :: bandfft_kpt_tab
 949  type(mpi_type),intent(inout) :: mpi_enreg
 950 !arrays
 951  integer,intent(in) :: ngfft(18)
 952  real(dp),intent(in) :: occ_k(nband_k)
 953  real(dp),intent(out) :: rhoaug(n4,n5,n6)
 954  real(dp),intent(in) :: cwavef(2,npw_k*blocksize)
 955  real(dp),target,intent(inout) :: wfraug(2,n4,n5,n6*ndat)
 956 
 957 !Local variables-------------------------------
 958 !scalars
 959  integer :: bandpp,bandpp_sym,ier,iibandpp,ikpt_this_proc,ind_occ,ind_occ1,ind_occ2,ipw
 960  integer :: istwf_k_,jjbandpp,me_fft,nd3,nproc_band,nproc_fft,npw_fft
 961  integer :: old_me_g0=0,spaceComm=0,tim_fourwf,use_gpu_cuda_
 962  integer,pointer :: idatarecv0,ndatarecv,ndatarecv_tot,ndatasend_sym
 963  logical :: flag_inv_sym,have_to_reequilibrate
 964  real(dp) :: weight,weight1,weight2
 965  type(bandfft_kpt_type),pointer :: bandfft_kpt_ptr
 966 !arrays
 967  integer,ABI_CONTIGUOUS pointer :: indices_pw_fft(:),kg_k_fft(:,:),kg_k_gather(:,:),kg_k_gather_sym(:,:)
 968  integer,ABI_CONTIGUOUS pointer :: rdispls(:),rdispls_sym(:)
 969  integer,ABI_CONTIGUOUS pointer :: recvcounts(:),recvcount_fft(:),recvcounts_sym(:),recvcounts_sym_tot(:)
 970  integer,ABI_CONTIGUOUS pointer :: recvdisp_fft(:),sdispls(:),sdispls_sym(:)
 971  integer,ABI_CONTIGUOUS pointer :: sendcounts(:),sendcount_fft(:),sendcounts_sym(:),sendcounts_sym_all(:)
 972  integer,ABI_CONTIGUOUS pointer :: senddisp_fft(:),tab_proc(:)
 973  integer,allocatable :: rdisplsloc(:)
 974  integer,allocatable :: recvcountsloc(:),sdisplsloc(:)
 975  integer,allocatable :: sendcountsloc(:)
 976  integer,allocatable :: index_wavef_band(:),index_wavef_send(:)
 977  integer,pointer :: gbound_(:,:)
 978  real(dp) :: dummy(2,1),tsec(2)
 979  real(dp),allocatable :: buff_wf(:,:),cwavef_alltoall1(:,:),cwavef_alltoall2(:,:)
 980  real(dp),allocatable :: cwavef_fft(:,:), cwavef_fft_tr(:,:)
 981  real(dp),allocatable :: weight_t(:),weight1_t(:),weight2_t(:)
 982  real(dp),pointer :: ewavef_alltoall_sym(:,:),wfraug_ptr(:,:,:,:)
 983 
 984 ! *************************************************************************
 985 
 986  ABI_CHECK((option_fourwf/=3),'Option=3 (FFT r->g) not implemented')
 987  ABI_CHECK((mpi_enreg%bandpp==ndat),'BUG: bandpp/=ndat')
 988 
 989  spaceComm=mpi_enreg%comm_band
 990  nproc_band = mpi_enreg%nproc_band
 991  nproc_fft  = mpi_enreg%nproc_fft
 992  bandpp     = mpi_enreg%bandpp
 993  me_fft     = mpi_enreg%me_fft
 994 
 995  use_gpu_cuda_=0;if (present(use_gpu_cuda)) use_gpu_cuda_=use_gpu_cuda
 996 
 997  if (present(bandfft_kpt_tab)) then
 998    bandfft_kpt_ptr => bandfft_kpt_tab
 999  else
1000    ikpt_this_proc=bandfft_kpt_get_ikpt()
1001    bandfft_kpt_ptr => bandfft_kpt(ikpt_this_proc)
1002  end if
1003 
1004  istwf_k_=istwf_k
1005  flag_inv_sym = (istwf_k_==2 .and. any(ngfft(7) == [401,402,312]))
1006  if (option_fourwf==0) flag_inv_sym=((flag_inv_sym).and.(use_gpu_cuda_==0))
1007 
1008  if (flag_inv_sym) then
1009    istwf_k_       = 1
1010    if (modulo(bandpp,2)==0) then
1011      bandpp_sym   = bandpp/2
1012    else
1013      bandpp_sym   = bandpp
1014    end if
1015  end if
1016 
1017 !====================================================================================
1018  ABI_ALLOCATE(sendcountsloc,(nproc_band))
1019  ABI_ALLOCATE(sdisplsloc   ,(nproc_band))
1020  ABI_ALLOCATE(recvcountsloc,(nproc_band))
1021  ABI_ALLOCATE(rdisplsloc   ,(nproc_band))
1022 
1023  recvcounts   =>bandfft_kpt_ptr%recvcounts(:)
1024  sendcounts   =>bandfft_kpt_ptr%sendcounts(:)
1025  rdispls      =>bandfft_kpt_ptr%rdispls   (:)
1026  sdispls      =>bandfft_kpt_ptr%sdispls   (:)
1027  ndatarecv    =>bandfft_kpt_ptr%ndatarecv
1028 
1029  kg_k_gather  =>bandfft_kpt_ptr%kg_k_gather(:,:)
1030  gbound_      =>bandfft_kpt_ptr%gbound(:,:)
1031 
1032  if (flag_inv_sym ) then
1033    idatarecv0           =>bandfft_kpt_ptr%idatarecv0
1034    ndatarecv_tot        =>bandfft_kpt_ptr%ndatarecv_tot
1035    ndatasend_sym        =>bandfft_kpt_ptr%ndatasend_sym
1036    kg_k_gather_sym      =>bandfft_kpt_ptr%kg_k_gather_sym(:,:)
1037    rdispls_sym          =>bandfft_kpt_ptr%rdispls_sym(:)
1038    recvcounts_sym       =>bandfft_kpt_ptr%recvcounts_sym(:)
1039    recvcounts_sym_tot   =>bandfft_kpt_ptr%recvcounts_sym_tot(:)
1040    sdispls_sym          =>bandfft_kpt_ptr%sdispls_sym(:)
1041    sendcounts_sym       =>bandfft_kpt_ptr%sendcounts_sym(:)
1042    sendcounts_sym_all   =>bandfft_kpt_ptr%sendcounts_sym_all(:)
1043    tab_proc             =>bandfft_kpt_ptr%tab_proc(:)
1044  end if
1045 
1046  ABI_ALLOCATE(cwavef_alltoall2,(2,ndatarecv*bandpp))
1047  if ( ((.not.flag_inv_sym) .and. (bandpp>1) ) .or. flag_inv_sym )then
1048    ABI_ALLOCATE(cwavef_alltoall1,(2,ndatarecv*bandpp))
1049  end if
1050 
1051  recvcountsloc(:)=recvcounts(:)*2*bandpp
1052  rdisplsloc(:)=rdispls(:)*2*bandpp
1053  sendcountsloc(:)=sendcounts(:)*2
1054  sdisplsloc(:)=sdispls(:)*2
1055 
1056  call timab(547,1,tsec)
1057  call xmpi_alltoallv(cwavef,sendcountsloc,sdisplsloc,cwavef_alltoall2,&
1058 & recvcountsloc,rdisplsloc,spaceComm,ier)
1059  call timab(547,2,tsec)
1060 
1061 !If me_fft==0, I have the G=0 vector, but keep for the record the old value
1062  if (me_fft==0) then
1063    old_me_g0=mpi_enreg%me_g0
1064    mpi_enreg%me_g0=1
1065  end if
1066 
1067  tim_fourwf=16
1068 
1069 !Eventually adjust load balancing for FFT (by changing FFT distrib)
1070  have_to_reequilibrate = bandfft_kpt_ptr%have_to_reequilibrate
1071  if(have_to_reequilibrate) then
1072    npw_fft =  bandfft_kpt_ptr%npw_fft
1073    sendcount_fft  => bandfft_kpt_ptr%sendcount_fft(:)
1074    recvcount_fft  => bandfft_kpt_ptr%recvcount_fft(:)
1075    senddisp_fft   => bandfft_kpt_ptr%senddisp_fft(:)
1076    recvdisp_fft   => bandfft_kpt_ptr%recvdisp_fft(:)
1077    indices_pw_fft => bandfft_kpt_ptr%indices_pw_fft(:)
1078    kg_k_fft       => bandfft_kpt_ptr%kg_k_fft(:,:)
1079    ABI_ALLOCATE( buff_wf, (2,ndatarecv*bandpp) ) ! for sorting cgwavef
1080    ABI_ALLOCATE( cwavef_fft, (2,npw_fft*bandpp) )
1081    if(bandpp>1) then
1082      ABI_ALLOCATE( cwavef_fft_tr, (2,npw_fft*bandpp))
1083    end if
1084  end if
1085 
1086  if (option_fourwf==0) wfraug(:,:,:,:)=zero
1087 
1088 !====================================================================
1089  if ((.not.(flag_inv_sym)) .and. (bandpp==1)) then
1090 
1091 !  Compute the index of the band
1092    ind_occ = (iblock-1)*blocksize + mpi_enreg%me_band + 1
1093 
1094    if(abs(occ_k(ind_occ))>=tol8.or.option_fourwf==0) then
1095 
1096 !    Compute the weight of the band
1097      weight=occ_k(ind_occ)*wtk/ucvol
1098 
1099      if(have_to_reequilibrate) then
1100 !      filling of sorted send buffers before exchange
1101        do ipw = 1 ,ndatarecv
1102          buff_wf(1:2, indices_pw_fft(ipw) ) = cwavef_alltoall2(1:2,ipw)
1103        end do
1104        call xmpi_alltoallv(buff_wf,2*sendcount_fft,2*senddisp_fft,  &
1105 &       cwavef_fft,2*recvcount_fft, 2*recvdisp_fft, mpi_enreg%comm_fft,ier)
1106        call fourwf(1,rhoaug,cwavef_fft,dummy,wfraug,gbound_,gbound_,&
1107 &       istwf_k_,kg_k_fft,kg_k_fft,mgfft,mpi_enreg,1,&
1108 &       ngfft,npw_fft,1,n4,n5,n6,option_fourwf,mpi_enreg%paral_kgb,tim_fourwf,weight,weight,&
1109 &       use_gpu_cuda=use_gpu_cuda_)
1110      else
1111        call fourwf(1,rhoaug,cwavef_alltoall2,dummy,wfraug,gbound_,gbound_,&
1112 &       istwf_k_,kg_k_gather,kg_k_gather,mgfft,mpi_enreg,1,&
1113 &       ngfft,ndatarecv,1,n4,n5,n6,option_fourwf,mpi_enreg%paral_kgb,tim_fourwf,weight,weight,&
1114 &       use_gpu_cuda=use_gpu_cuda_)
1115      end if
1116      if (option_fourwf==0.and.nproc_fft>1) then
1117        if (me_fft>0) then
1118          nd3=(ngfft(3)-1)/nproc_fft+1
1119          wfraug(:,:,:,me_fft*nd3+1:me_fft*nd3+nd3)=wfraug(:,:,:,1:nd3)
1120          wfraug(:,:,:,1:nd3)=zero
1121        end if
1122        call xmpi_sum(wfraug,mpi_enreg%comm_fft,ier)
1123      end if
1124    end if
1125 
1126 !====================================================================
1127  else if ((.not.(flag_inv_sym)) .and. (bandpp>1) ) then
1128 
1129 !  -------------------------------------------------------------
1130 !  Computation of the index to class the waves functions below bandpp
1131 !  -------------------------------------------------------------
1132    call prep_index_wavef_bandpp(nproc_band,bandpp,&
1133 &   1,ndatarecv,&
1134 &   recvcounts,rdispls,&
1135 &   index_wavef_band)
1136 
1137 !  -------------------------------------------------------
1138 !  Sorting of the wave functions below bandpp
1139 !  -------------------------------------------------------
1140    cwavef_alltoall1(:,:) = cwavef_alltoall2(:,index_wavef_band)
1141 
1142    if(have_to_reequilibrate) then
1143 !    filling of sorted send buffers before exchange
1144      do iibandpp=1,bandpp
1145        do ipw = 1 ,ndatarecv
1146          buff_wf(1:2, iibandpp + bandpp*(indices_pw_fft(ipw)-1)) = &
1147 &         cwavef_alltoall1(1:2,ipw + ndatarecv*(iibandpp-1))
1148        end do
1149      end do
1150      call xmpi_alltoallv(buff_wf,2*bandpp*sendcount_fft,2*bandpp*senddisp_fft,  &
1151 &     cwavef_fft_tr,2*bandpp*recvcount_fft, 2*bandpp*recvdisp_fft, mpi_enreg%comm_fft,ier)
1152      do iibandpp=1,bandpp
1153        do ipw = 1 ,npw_fft
1154          cwavef_fft(1:2,  ipw + npw_fft*(iibandpp -1)) = cwavef_fft_tr(1:2,  iibandpp + bandpp*(ipw-1))
1155        end do
1156      end do
1157    end if
1158 
1159 !  -------------------
1160 !  Fourier calculation
1161 !  -------------------
1162 !  Cuda version
1163    if(use_gpu_cuda_==1) then
1164      ABI_ALLOCATE(weight_t,(bandpp))
1165      do iibandpp=1,bandpp
1166 !      Compute the index of the band
1167        ind_occ = (iblock-1)*blocksize + (mpi_enreg%me_band * bandpp) + iibandpp
1168 !      Compute the weight of the band
1169        weight_t(iibandpp)=occ_k(ind_occ)*wtk/ucvol
1170        if(abs(occ_k(ind_occ)) < tol8) weight_t(iibandpp) = zero
1171      end do
1172 !    Accumulate time because it is not done in gpu_fourwf
1173      call timab(240+tim_fourwf,1,tsec)
1174 #if defined HAVE_GPU_CUDA
1175      call gpu_fourwf(1,rhoaug,&
1176 &     cwavef_alltoall1,&
1177 &     dummy,wfraug,gbound_,gbound_,&
1178 &     istwf_k_,kg_k_gather,kg_k_gather,mgfft,mpi_enreg,bandpp,&
1179 &     ngfft,ndatarecv,1,n4,n5,n6,option_fourwf,mpi_enreg%paral_kgb,&
1180 &     tim_fourwf,weight_t,weight_t)
1181 #endif
1182      call timab(240+tim_fourwf,2,tsec)
1183      ABI_DEALLOCATE(weight_t)
1184 
1185 !  Standard version
1186    else
1187      do iibandpp=1,bandpp
1188 !      Compute the index of the band
1189        ind_occ = (iblock-1)*blocksize + (mpi_enreg%me_band * bandpp) + iibandpp
1190 !      Compute the weight of the band
1191        weight=occ_k(ind_occ)*wtk/ucvol
1192        if (option_fourwf==0) then
1193          wfraug_ptr => wfraug(:,:,:,(iibandpp-1)*n6+1:iibandpp*n6)
1194        else
1195          wfraug_ptr => wfraug
1196        end if
1197        if (abs(occ_k(ind_occ)) >=tol8.or.option_fourwf==0) then
1198          if(have_to_reequilibrate) then
1199            call fourwf(1,rhoaug, &
1200 &           cwavef_fft(:,(npw_fft*(iibandpp-1))+1:(npw_fft*iibandpp)), &
1201 &           dummy,wfraug_ptr,gbound_,gbound_,&
1202 &           istwf_k_,kg_k_fft,kg_k_fft,mgfft,mpi_enreg,1,&
1203 &           ngfft,npw_fft,1,n4,n5,n6,option_fourwf,mpi_enreg%paral_kgb,tim_fourwf,weight,weight,&
1204 &           use_gpu_cuda=use_gpu_cuda_)
1205          else
1206            call fourwf(1,rhoaug,&
1207 &           cwavef_alltoall1(:,(ndatarecv*(iibandpp-1))+1:(ndatarecv*iibandpp)),&
1208 &           dummy,wfraug_ptr,gbound_,gbound_,&
1209 &           istwf_k_,kg_k_gather,kg_k_gather,mgfft,mpi_enreg,1,&
1210 &           ngfft,ndatarecv,1,n4,n5,n6,option_fourwf,mpi_enreg%paral_kgb,&
1211 &           tim_fourwf,weight,weight)
1212          end if
1213          if (option_fourwf==0.and.nproc_fft>1) then
1214            if (me_fft>0) then
1215              nd3=(ngfft(3)-1)/nproc_fft+1
1216              wfraug_ptr(:,:,:,me_fft*nd3+1:me_fft*nd3+nd3)=wfraug_ptr(:,:,:,1:nd3)
1217              wfraug_ptr(:,:,:,1:nd3)=zero
1218            end if
1219            call xmpi_sum(wfraug_ptr,mpi_enreg%comm_fft,ier)
1220          end if
1221        end if
1222      end do
1223    end if ! (use_gpu_cuda==1)
1224 
1225 !  -----------------------------------------------------
1226 !  Sorting waves functions below the processors
1227 !  -----------------------------------------------------
1228 !  cwavef_alltoall(:,index_wavef_band) = cwavef_alltoall(:,:)   ! NOT NEEDED
1229    ABI_DEALLOCATE(index_wavef_band)
1230 
1231 !====================================================================
1232  else if (flag_inv_sym) then
1233 
1234 !  -------------------------------------------------------------
1235 !  Computation of the index to class the waves functions below bandpp
1236 !  -------------------------------------------------------------
1237    call prep_index_wavef_bandpp(nproc_band,bandpp,&
1238 &   1,ndatarecv,&
1239 &   recvcounts,rdispls,&
1240 &   index_wavef_band)
1241 
1242 !  -------------------------------------------------------
1243 !  Sorting the wave functions below bandpp
1244 !  -------------------------------------------------------
1245    cwavef_alltoall1(:,:) = cwavef_alltoall2(:,index_wavef_band)
1246 
1247 !  ------------------------------------------------------------
1248 !  We associate the waves functions by two
1249 !  ------------------------------------------------------------
1250    call prep_wavef_sym_do(mpi_enreg,bandpp,1,&
1251 &   ndatarecv,&
1252 &   ndatarecv_tot,ndatasend_sym,tab_proc,&
1253 &   cwavef_alltoall1,&
1254 &   sendcounts_sym,sdispls_sym,&
1255 &   recvcounts_sym,rdispls_sym,&
1256 &   ewavef_alltoall_sym,&
1257 &   index_wavef_send)
1258 
1259 !  ------------------------------------------------------------
1260 !  Fourier calculation
1261 !  ------------------------------------------------------------
1262 !  Cuda version
1263    if (use_gpu_cuda_==1) then
1264      ABI_ALLOCATE(weight1_t,(bandpp_sym))
1265      ABI_ALLOCATE(weight2_t,(bandpp_sym))
1266      do iibandpp=1,bandpp_sym
1267        if (bandpp/=1) then
1268          ind_occ1 = (iblock-1)*blocksize + (mpi_enreg%me_band * bandpp) + (2*iibandpp-1)
1269          ind_occ2 = (iblock-1)*blocksize + (mpi_enreg%me_band * bandpp) + (2*iibandpp  )
1270        else
1271          ind_occ1 = (iblock-1)*blocksize + (mpi_enreg%me_band * bandpp) + 1
1272          ind_occ2 = ind_occ1
1273        end if
1274        weight1_t(iibandpp) = occ_k(ind_occ1)*wtk/ucvol
1275        weight2_t(iibandpp) = occ_k(ind_occ2)*wtk/ucvol
1276      end do
1277      call timab(240+tim_fourwf,1,tsec)
1278 #if defined HAVE_GPU_CUDA
1279      call gpu_fourwf(1,rhoaug,&
1280 &     ewavef_alltoall_sym,&
1281 &     dummy,wfraug,gbound_,gbound_,&
1282 &     istwf_k_,kg_k_gather_sym,kg_k_gather_sym,mgfft,mpi_enreg,bandpp_sym,&
1283 &     ngfft,ndatarecv_tot,1,n4,n5,n6,option_fourwf,mpi_enreg%paral_kgb,&
1284 &     tim_fourwf,weight1_t,weight2_t)
1285 #endif
1286      call timab(240+tim_fourwf,2,tsec)
1287      ABI_DEALLOCATE(weight1_t)
1288      ABI_DEALLOCATE(weight2_t)
1289 
1290 !  Standard version
1291    else
1292      if (option_fourwf==0.and.bandpp>1) then
1293        ABI_ALLOCATE(wfraug_ptr,(2,n4,n5,n6))
1294      else
1295        wfraug_ptr => wfraug
1296      end if
1297      do iibandpp=1,bandpp_sym
1298        if (bandpp/=1) then
1299          ind_occ1 = (iblock-1)*blocksize + (mpi_enreg%me_band * bandpp) + (2*iibandpp-1)
1300          ind_occ2 = (iblock-1)*blocksize + (mpi_enreg%me_band * bandpp) + (2*iibandpp  )
1301        else
1302          ind_occ1 = (iblock-1)*blocksize + (mpi_enreg%me_band * bandpp) + 1
1303          ind_occ2 = ind_occ1
1304        end if
1305        weight1 = occ_k(ind_occ1)*wtk/ucvol
1306        weight2 = occ_k(ind_occ2)*wtk/ucvol
1307        call fourwf(1,rhoaug,&
1308 &       ewavef_alltoall_sym(:,(ndatarecv_tot*(iibandpp-1))+1:(ndatarecv_tot*iibandpp)),&
1309 &       dummy,wfraug_ptr,gbound_,gbound_,&
1310 &       istwf_k_,kg_k_gather_sym,kg_k_gather_sym,mgfft,mpi_enreg,1,&
1311 &       ngfft,ndatarecv_tot,1,n4,n5,n6,option_fourwf,mpi_enreg%paral_kgb,&
1312 &       tim_fourwf,weight1,weight2)
1313        if (option_fourwf==0) then
1314          if (modulo(bandpp,2)==0) then
1315            jjbandpp=2*iibandpp-1
1316            wfraug(1,:,:,(jjbandpp-1)*n6+1:jjbandpp*n6)=wfraug_ptr(1,:,:,1:n6)
1317            wfraug(1,:,:,(jjbandpp)*n6+1:(jjbandpp+1)*n6)=wfraug_ptr(2,:,:,1:n6)
1318          else if (bandpp>1) then
1319            wfraug(1,:,:,(iibandpp-1)*n6+1:iibandpp*n6)=wfraug_ptr(1,:,:,1:n6)
1320          end if
1321          if (nproc_fft>1) then
1322            if (me_fft>0) then
1323              nd3=(ngfft(3)-1)/nproc_fft+1
1324              wfraug(1,:,:,me_fft*nd3+1:me_fft*nd3+nd3)=wfraug(1,:,:,1:nd3)
1325              wfraug(1,:,:,1:nd3)=zero
1326            end if
1327            call xmpi_sum(wfraug,mpi_enreg%comm_fft,ier)
1328          end if
1329        end if
1330      end do
1331      if (option_fourwf==0.and.bandpp>1) then
1332        ABI_DEALLOCATE(wfraug_ptr)
1333      end if
1334    end if ! (use_gpu_cuda==1)
1335 
1336 !  ------------------------------------------------------------
1337 !  We dissociate each wave function in two waves functions
1338 !  gwavef is classed below of bandpp
1339 !  ------------------------------------------------------------
1340    call prep_wavef_sym_undo(mpi_enreg,bandpp,1,&
1341 &   ndatarecv,&
1342 &   ndatarecv_tot,ndatasend_sym,idatarecv0,&
1343 &   cwavef_alltoall1,&
1344 &   sendcounts_sym,sdispls_sym,&
1345 &   recvcounts_sym,rdispls_sym,&
1346 &   ewavef_alltoall_sym,&
1347 &   index_wavef_send)
1348 
1349    ABI_DEALLOCATE(ewavef_alltoall_sym)
1350    ABI_DEALLOCATE(index_wavef_send)
1351 
1352 !  -------------------------------------------------------
1353 !  Sorting waves functions below the processors
1354 !  -------------------------------------------------------
1355 !  cwavef_alltoall(:,index_wavef_band) = cwavef_alltoall(:,:) ! NOT NEEDED
1356 
1357    ABI_DEALLOCATE(index_wavef_band)
1358 
1359  end if
1360 
1361 !====================================================================
1362  if (me_fft==0) mpi_enreg%me_g0=old_me_g0
1363  if(have_to_reequilibrate) then
1364    ABI_DEALLOCATE(buff_wf)
1365    ABI_DEALLOCATE(cwavef_fft)
1366    if(bandpp > 1) then
1367      ABI_DEALLOCATE(cwavef_fft_tr)
1368    end if
1369  end if
1370  ABI_DEALLOCATE(sendcountsloc)
1371  ABI_DEALLOCATE(sdisplsloc)
1372  ABI_DEALLOCATE(recvcountsloc)
1373  ABI_DEALLOCATE(rdisplsloc)
1374  ABI_DEALLOCATE(cwavef_alltoall2)
1375  if ( ((.not.flag_inv_sym) .and. (bandpp>1) ) .or. flag_inv_sym ) then
1376    ABI_DEALLOCATE(cwavef_alltoall1)
1377  end if
1378 
1379 end subroutine prep_fourwf

ABINIT/prep_getghc [ Functions ]

[ Top ] [ Functions ]

NAME

 prep_getghc

FUNCTION

 this routine prepares the data to the call of getghc.

INPUTS

  blocksize= size of block for FFT
  cpopt=flag defining the status of cprjin%cp(:)=<Proj_i|Cnk> scalars (see below, side effects)
  cwavef(2,npw*my_nspinor*blocksize)=planewave coefficients of wavefunction.
  gs_hamk <type(gs_hamiltonian_type)>=all data for the hamiltonian at k
  gvnlc=matrix elements <G|Vnonlocal|C>
  lambda=factor to be used when computing <G|H-lambda.S|C> - only for sij_opt=-1
         Typically lambda is the eigenvalue (or its guess)
  mpi_enreg=informations about mpi parallelization
  prtvol=control print volume and debugging output
  sij_opt= -PAW ONLY-  if  0, only matrix elements <G|H|C> have to be computed
     (S=overlap)       if  1, matrix elements <G|S|C> have to be computed in gsc in addition to ghc
                       if -1, matrix elements <G|H-lambda.S|C> have to be computed in ghc (gsc not used)

OUTPUT

  gwavef=(2,npw*my_nspinor*blocksize)=matrix elements <G|H|C> (if sij_opt>=0)
                                  or <G|H-lambda.S|C> (if sij_opt=-1).
  swavef=(2,npw*my_nspinor*blocksize)=matrix elements <G|S|C>.

SIDE EFFECTS

  ====== if gs_hamk%usepaw==1
  cwaveprj(natom,my_nspinor*bandpp)= wave functions at k projected with nl projectors

PARENTS

      chebfi,lobpcgwf,m_lobpcgwf,mkresi

CHILDREN

      dcopy,multithreaded_getghc,prep_index_wavef_bandpp,prep_sort_wavef_spin
      prep_wavef_sym_do,prep_wavef_sym_undo,timab,xmpi_alltoallv

SOURCE

 99 subroutine prep_getghc(cwavef,gs_hamk,gvnlc,gwavef,swavef,lambda,blocksize,&
100 &                      mpi_enreg,prtvol,sij_opt,cpopt,cwaveprj,&
101 &                      already_transposed) ! optional argument
102 
103 
104 !This section has been created automatically by the script Abilint (TD).
105 !Do not modify the following lines by hand.
106 #undef ABI_FUNC
107 #define ABI_FUNC 'prep_getghc'
108 !End of the abilint section
109 
110  implicit none
111 
112 !Arguments ------------------------------------
113 !scalars
114  integer,intent(in) :: blocksize,cpopt,prtvol,sij_opt
115  logical, intent(in),optional :: already_transposed
116  real(dp),intent(in) :: lambda
117  type(gs_hamiltonian_type),intent(inout) :: gs_hamk
118  type(mpi_type),intent(inout) :: mpi_enreg
119 !arrays
120  real(dp),intent(inout) :: cwavef(:,:),gvnlc (:,:),gwavef(:,:),swavef(:,:)
121  type(pawcprj_type), intent(inout) :: cwaveprj(:,:)
122 
123 !Local variables-------------------------------
124 !scalars
125  integer,parameter :: tim_getghc=6
126  integer :: bandpp,bandpp_sym,idatarecv0,ier,ikpt_this_proc,iscalc,mcg,my_nspinor
127  integer :: nbval,ndatarecv,ndatarecv_tot,ndatasend_sym,nproc_band,nproc_fft
128  integer :: old_me_g0,spaceComm=0
129  logical :: flag_inv_sym, do_transpose
130  character(len=100) :: msg
131 !arrays
132  integer,allocatable :: index_wavef_band(:),index_wavef_send(:),index_wavef_spband(:)
133  integer,allocatable :: rdisplsloc(:),recvcountsloc(:),sdisplsloc(:),sendcountsloc(:)
134  integer,ABI_CONTIGUOUS pointer :: kg_k_gather_sym(:,:)
135  integer,ABI_CONTIGUOUS pointer :: rdispls(:),rdispls_sym(:)
136  integer,ABI_CONTIGUOUS pointer :: recvcounts(:),recvcounts_sym(:),recvcounts_sym_tot(:)
137  integer,ABI_CONTIGUOUS pointer :: sdispls(:),sdispls_sym(:)
138  integer,ABI_CONTIGUOUS pointer :: sendcounts(:),sendcounts_sym(:),sendcounts_sym_all(:)
139  integer,ABI_CONTIGUOUS pointer :: tab_proc(:)
140  real(dp) :: tsec(2)
141  real(dp),allocatable,target :: cwavef_alltoall1(:,:),gvnlc_alltoall1(:,:)
142  real(dp),allocatable,target :: gwavef_alltoall1(:,:),swavef_alltoall1(:,:)
143  real(dp),allocatable,target :: cwavef_alltoall2(:,:),gvnlc_alltoall2(:,:)
144  real(dp),allocatable,target :: gwavef_alltoall2(:,:),swavef_alltoall2(:,:)
145  real(dp),pointer :: ewavef_alltoall_sym(:,:),gvnlc_alltoall_sym(:,:)
146  real(dp),pointer :: gwavef_alltoall_sym(:,:)
147  real(dp),pointer :: swavef_alltoall_sym(:,:)
148 
149 ! *************************************************************************
150 
151  call timab(630,1,tsec)
152  call timab(631,3,tsec)
153 
154 !Some inits
155  nproc_band = mpi_enreg%nproc_band
156  nproc_fft  = mpi_enreg%nproc_fft
157  bandpp     = mpi_enreg%bandpp
158  my_nspinor = max(1,gs_hamk%nspinor/mpi_enreg%nproc_spinor)
159 
160  do_transpose = .true.
161  if(present(already_transposed)) then
162    if(already_transposed) do_transpose = .false.
163  end if
164 
165  flag_inv_sym = (gs_hamk%istwf_k==2 .and. any(gs_hamk%ngfft(7) == [401,402,312]))
166  if (flag_inv_sym) then
167    gs_hamk%istwf_k = 1
168    if (modulo(bandpp,2)==0) bandpp_sym = bandpp/2
169    if (modulo(bandpp,2)/=0) bandpp_sym = bandpp
170  end if
171 
172 !Check sizes
173  mcg=2*gs_hamk%npw_fft_k*my_nspinor*bandpp
174  if (do_transpose) mcg=2*gs_hamk%npw_k*my_nspinor*blocksize
175  if (size(cwavef)<mcg) then
176    msg='wrong size for cwavef!'
177    MSG_BUG(msg)
178  end if
179  if (size(gwavef)<mcg) then
180    msg='wrong size for gwavef!'
181    MSG_BUG(msg)
182  end if
183  if (size(gvnlc)<mcg) then
184    msg='wrong size for gvnlc!'
185    MSG_BUG(msg)
186  end if
187  if (sij_opt==1) then
188    if (size(swavef)<mcg) then
189      msg='wrong size for swavef!'
190      MSG_BUG(msg)
191    end if
192  end if
193  if (gs_hamk%usepaw==1.and.cpopt>=0) then
194    if (size(cwaveprj)<gs_hamk%natom*my_nspinor*bandpp) then
195      msg='wrong size for cwaveprj!'
196      MSG_BUG(msg)
197    end if
198  end if
199 
200 !====================================================================================
201 
202  spaceComm=mpi_enreg%comm_fft
203  if(mpi_enreg%paral_kgb==1) spaceComm=mpi_enreg%comm_band
204 
205  ikpt_this_proc=bandfft_kpt_get_ikpt()
206 
207  ABI_ALLOCATE(sendcountsloc,(nproc_band))
208  ABI_ALLOCATE(sdisplsloc   ,(nproc_band))
209  ABI_ALLOCATE(recvcountsloc,(nproc_band))
210  ABI_ALLOCATE(rdisplsloc   ,(nproc_band))
211 
212  recvcounts   =>bandfft_kpt(ikpt_this_proc)%recvcounts(:)
213  sendcounts   =>bandfft_kpt(ikpt_this_proc)%sendcounts(:)
214  rdispls      =>bandfft_kpt(ikpt_this_proc)%rdispls   (:)
215  sdispls      =>bandfft_kpt(ikpt_this_proc)%sdispls   (:)
216  ndatarecv    = bandfft_kpt(ikpt_this_proc)%ndatarecv
217 
218  if (flag_inv_sym ) then
219    idatarecv0           = bandfft_kpt(ikpt_this_proc)%idatarecv0
220    ndatarecv_tot        = bandfft_kpt(ikpt_this_proc)%ndatarecv_tot
221    ndatasend_sym        = bandfft_kpt(ikpt_this_proc)%ndatasend_sym
222    kg_k_gather_sym      =>bandfft_kpt(ikpt_this_proc)%kg_k_gather_sym(:,:)
223    rdispls_sym          =>bandfft_kpt(ikpt_this_proc)%rdispls_sym(:)
224    recvcounts_sym       =>bandfft_kpt(ikpt_this_proc)%recvcounts_sym(:)
225    recvcounts_sym_tot   =>bandfft_kpt(ikpt_this_proc)%recvcounts_sym_tot(:)
226    sdispls_sym          =>bandfft_kpt(ikpt_this_proc)%sdispls_sym(:)
227    sendcounts_sym       =>bandfft_kpt(ikpt_this_proc)%sendcounts_sym(:)
228    sendcounts_sym_all   =>bandfft_kpt(ikpt_this_proc)%sendcounts_sym_all(:)
229    tab_proc             =>bandfft_kpt(ikpt_this_proc)%tab_proc(:)
230  end if
231  iscalc=(sij_opt+1)/2  ! 0 if S not calculated, 1 otherwise
232  nbval=(ndatarecv*my_nspinor*bandpp)*iscalc
233 
234  if ( ((.not.flag_inv_sym) .and. bandpp==1 .and. mpi_enreg%paral_spinor==0 .and. my_nspinor==2 ).or. &
235 & ((.not.flag_inv_sym) .and. bandpp>1) .or.  flag_inv_sym  ) then
236    ABI_ALLOCATE(cwavef_alltoall1,(2,ndatarecv*my_nspinor*bandpp))
237    ABI_ALLOCATE(gwavef_alltoall1,(2,ndatarecv*my_nspinor*bandpp))
238    ABI_ALLOCATE(swavef_alltoall1,(2,ndatarecv*my_nspinor*bandpp))
239    ABI_ALLOCATE(gvnlc_alltoall1,(2,ndatarecv*my_nspinor*bandpp))
240    swavef_alltoall1(:,:)=zero
241    gvnlc_alltoall1(:,:)=zero
242    cwavef_alltoall1(:,:)=zero
243    gwavef_alltoall1(:,:)=zero
244  end if
245  ABI_ALLOCATE(cwavef_alltoall2,(2,ndatarecv*my_nspinor*bandpp))
246  ABI_ALLOCATE(gwavef_alltoall2,(2,ndatarecv*my_nspinor*bandpp))
247  ABI_ALLOCATE(swavef_alltoall2,(2,ndatarecv*my_nspinor*bandpp))
248  ABI_ALLOCATE(gvnlc_alltoall2,(2,ndatarecv*my_nspinor*bandpp))
249  swavef_alltoall2(:,:)=zero
250  gvnlc_alltoall2(:,:)=zero
251  cwavef_alltoall2(:,:)=zero
252  gwavef_alltoall2(:,:)=zero
253 
254  recvcountsloc(:)=recvcounts(:)*2*my_nspinor*bandpp
255  rdisplsloc(:)=rdispls(:)*2*my_nspinor*bandpp
256  sendcountsloc(:)=sendcounts(:)*2*my_nspinor
257  sdisplsloc(:)=sdispls(:)*2*my_nspinor
258  call timab(631,2,tsec)
259 
260  if(do_transpose) then
261    call timab(545,3,tsec)
262    if ( ((.not.flag_inv_sym) .and. bandpp==1 .and. mpi_enreg%paral_spinor==0 .and. my_nspinor==2 ).or. &
263 &   ((.not.flag_inv_sym) .and. bandpp>1) .or.  flag_inv_sym  ) then
264      call xmpi_alltoallv(cwavef,sendcountsloc,sdisplsloc,cwavef_alltoall1,&
265 &     recvcountsloc,rdisplsloc,spaceComm,ier)
266    else
267      call xmpi_alltoallv(cwavef,sendcountsloc,sdisplsloc,cwavef_alltoall2,&
268 &     recvcountsloc,rdisplsloc,spaceComm,ier)
269    end if
270    call timab(545,2,tsec)
271  else
272    ! Here, we cheat, and use DCOPY to bypass some compiler's overzealous bound-checking
273    ! (ndatarecv*my_nspinor*bandpp might be greater than the declared size of cwavef)
274    call DCOPY(2*ndatarecv*my_nspinor*bandpp, cwavef, 1, cwavef_alltoall2, 1)
275  end if
276 
277  if(gs_hamk%istwf_k==2) then
278    old_me_g0=mpi_enreg%me_g0
279    if (mpi_enreg%me_fft==0) then
280      mpi_enreg%me_g0=1
281    else
282      mpi_enreg%me_g0=0
283    end if
284  end if
285 
286 !====================================================================
287  if ((.not.(flag_inv_sym)) .and. (bandpp==1)) then
288    if (do_transpose .and. mpi_enreg%paral_spinor==0.and.my_nspinor==2)then
289      call timab(632,3,tsec)
290 !    Sort to have all ispinor=1 first, then all ispinor=2
291      call prep_sort_wavef_spin(nproc_band,my_nspinor,ndatarecv,recvcounts,rdispls,index_wavef_spband)
292      cwavef_alltoall2(:,:)=cwavef_alltoall1(:,index_wavef_spband)
293      call timab(632,2,tsec)
294    end if
295 
296    call timab(635,3,tsec)
297    call multithreaded_getghc(cpopt,cwavef_alltoall2,cwaveprj,gwavef_alltoall2,swavef_alltoall2(:,1:nbval),&
298 &   gs_hamk,gvnlc_alltoall2,lambda,mpi_enreg,1,prtvol,sij_opt,tim_getghc,0)
299    call timab(635,2,tsec)
300 
301    if (do_transpose .and. mpi_enreg%paral_spinor==0.and.my_nspinor==2)then
302      call timab(634,3,tsec)
303      gwavef_alltoall1(:,index_wavef_spband)=gwavef_alltoall2(:,:)
304      if (sij_opt==1) swavef_alltoall1(:,index_wavef_spband)=swavef_alltoall2(:,:)
305      gvnlc_alltoall1(:,index_wavef_spband)=gvnlc_alltoall2(:,:)
306      ABI_DEALLOCATE(index_wavef_spband)
307      call timab(634,2,tsec)
308    end if
309 
310  else if ((.not.(flag_inv_sym)) .and. (bandpp>1)) then
311 !  -------------------------------------------------------------
312 !  Computation of the index to class the waves functions below bandpp
313 !  -------------------------------------------------------------
314 
315    if(do_transpose) then
316      call timab(632,3,tsec)
317      call prep_index_wavef_bandpp(nproc_band,bandpp,&
318 &     my_nspinor,ndatarecv, recvcounts,rdispls, index_wavef_band)
319 !  -------------------------------------------------------
320 !  Sorting of the waves functions below bandpp
321 !  -------------------------------------------------------
322      cwavef_alltoall2(:,:) = cwavef_alltoall1(:,index_wavef_band)
323      call timab(632,2,tsec)
324    end if
325 
326 !  ----------------------
327 !  Fourier transformation
328 !  ----------------------
329    call timab(636,3,tsec)
330    call multithreaded_getghc(cpopt,cwavef_alltoall2,cwaveprj,gwavef_alltoall2,swavef_alltoall2,gs_hamk,&
331 &   gvnlc_alltoall2,lambda,mpi_enreg,bandpp,prtvol,sij_opt,tim_getghc,0)
332    call timab(636,2,tsec)
333 
334 !  -----------------------------------------------------
335 !  Sorting of waves functions below the processors
336 !  -----------------------------------------------------
337    if(do_transpose) then
338      call timab(634,3,tsec)
339      gwavef_alltoall1(:,index_wavef_band) = gwavef_alltoall2(:,:)
340      if (sij_opt==1) swavef_alltoall1(:,index_wavef_band) = swavef_alltoall2(:,:)
341      gvnlc_alltoall1(:,index_wavef_band)  = gvnlc_alltoall2(:,:)
342      ABI_DEALLOCATE(index_wavef_band)
343      call timab(634,2,tsec)
344    end if
345 
346 
347  else if (flag_inv_sym) then
348 
349 !  -------------------------------------------------------------
350 !  Computation of the index to class the waves functions below bandpp
351 !  -------------------------------------------------------------
352    if(do_transpose) then
353      call timab(632,3,tsec)
354      call prep_index_wavef_bandpp(nproc_band,bandpp,&
355 &     my_nspinor,ndatarecv,&
356 &     recvcounts,rdispls,&
357 &     index_wavef_band)
358 
359 !  -------------------------------------------------------
360 !  Sorting the wave functions below bandpp
361 !  -------------------------------------------------------
362      cwavef_alltoall2(:,:) = cwavef_alltoall1(:,index_wavef_band)
363    end if
364 
365 !  ------------------------------------------------------------
366 !  We associate the waves functions by two
367 !  ------------------------------------------------------------
368    call prep_wavef_sym_do(mpi_enreg,bandpp,my_nspinor,&
369 &   ndatarecv,&
370 &   ndatarecv_tot,ndatasend_sym,tab_proc,&
371 &   cwavef_alltoall2,&
372 &   sendcounts_sym,sdispls_sym,&
373 &   recvcounts_sym,rdispls_sym,&
374 &   ewavef_alltoall_sym,&
375 &   index_wavef_send)
376 
377 !  ------------------------------------------------------------
378 !  Allocation
379 !  ------------------------------------------------------------
380    ABI_ALLOCATE(gwavef_alltoall_sym,(2,ndatarecv_tot*bandpp_sym))
381    ABI_ALLOCATE(swavef_alltoall_sym,(2,(ndatarecv_tot*bandpp_sym)*iscalc))
382    ABI_ALLOCATE(gvnlc_alltoall_sym ,(2,ndatarecv_tot*bandpp_sym))
383 
384    gwavef_alltoall_sym(:,:)=zero
385    swavef_alltoall_sym(:,:)=zero
386    gvnlc_alltoall_sym(:,:)=zero
387 
388    call timab(632,2,tsec)
389 
390 !  ------------------------------------------------------------
391 !  Fourier calculation
392 !  ------------------------------------------------------------
393    call timab(637,3,tsec)
394    call multithreaded_getghc(cpopt,ewavef_alltoall_sym,cwaveprj,gwavef_alltoall_sym,swavef_alltoall_sym,gs_hamk,&
395 &   gvnlc_alltoall_sym,lambda,mpi_enreg,bandpp_sym,prtvol,sij_opt,tim_getghc,1,&
396 &   kg_fft_k=kg_k_gather_sym)
397    call timab(637,2,tsec)
398 
399    call timab(633,3,tsec)
400 
401 !  ------------------------------------------------------------
402 !  We dissociate each wave function in two waves functions
403 !  gwavef is classed below of bandpp
404 !  ------------------------------------------------------------
405    call prep_wavef_sym_undo(mpi_enreg,bandpp,my_nspinor,&
406 &   ndatarecv,&
407 &   ndatarecv_tot,ndatasend_sym,idatarecv0,&
408 &   gwavef_alltoall2,&
409 &   sendcounts_sym,sdispls_sym,&
410 &   recvcounts_sym,rdispls_sym,&
411 &   gwavef_alltoall_sym,&
412 &   index_wavef_send)
413    if (sij_opt==1)then
414      call prep_wavef_sym_undo(mpi_enreg,bandpp,my_nspinor,&
415 &     ndatarecv,&
416 &     ndatarecv_tot,ndatasend_sym,idatarecv0,&
417 &     swavef_alltoall2,&
418 &     sendcounts_sym,sdispls_sym,&
419 &     recvcounts_sym,rdispls_sym,&
420 &     swavef_alltoall_sym,&
421 &     index_wavef_send)
422    end if
423    call prep_wavef_sym_undo(mpi_enreg,bandpp,my_nspinor,&
424 &   ndatarecv,&
425 &   ndatarecv_tot,ndatasend_sym,idatarecv0,&
426 &   gvnlc_alltoall2,&
427 &   sendcounts_sym,sdispls_sym,&
428 &   recvcounts_sym,rdispls_sym,&
429 &   gvnlc_alltoall_sym,&
430 &   index_wavef_send)
431 
432    ABI_DEALLOCATE(ewavef_alltoall_sym)
433    ABI_DEALLOCATE(index_wavef_send)
434    ABI_DEALLOCATE(gwavef_alltoall_sym)
435    ABI_DEALLOCATE(swavef_alltoall_sym)
436    ABI_DEALLOCATE(gvnlc_alltoall_sym)
437 
438 !  -------------------------------------------
439 !  We call getghc to calculate the nl matrix elements.
440 !  --------------------------------------------
441    gs_hamk%istwf_k=2
442    !!write(std_out,*)"Setting iswfk_k to 2"
443 
444    old_me_g0=mpi_enreg%me_g0
445    if (mpi_enreg%me_fft==0) then
446      mpi_enreg%me_g0=1
447    else
448      mpi_enreg%me_g0=0
449    end if
450    call timab(633,2,tsec)
451 
452    call timab(638,3,tsec)
453    call multithreaded_getghc(cpopt,cwavef_alltoall2,cwaveprj,gwavef_alltoall2,swavef_alltoall2,gs_hamk,&
454 &   gvnlc_alltoall2,lambda,mpi_enreg,bandpp,prtvol,sij_opt,tim_getghc,2)
455    call timab(638,2,tsec)
456 
457    call timab(634,3,tsec)
458    mpi_enreg%me_g0=old_me_g0
459 
460    gs_hamk%istwf_k=1
461 
462 !  -------------------------------------------------------
463 !  Sorting the wave functions below the processors
464 !  -------------------------------------------------------
465    if(do_transpose) then
466 !    cwavef_alltoall(:,index_wavef_band) = cwavef_alltoall(:,:)   ! NOT NEEDED
467      gwavef_alltoall1(:,index_wavef_band) = gwavef_alltoall2(:,:)
468      if (sij_opt==1) swavef_alltoall1(:,index_wavef_band) = swavef_alltoall2(:,:)
469      gvnlc_alltoall1(:,index_wavef_band)  = gvnlc_alltoall2(:,:)
470      ABI_DEALLOCATE(index_wavef_band)
471      call timab(634,2,tsec)
472    end if
473 
474  end if
475 !====================================================================
476 
477  if (gs_hamk%istwf_k==2) mpi_enreg%me_g0=old_me_g0
478 
479  if(do_transpose) then
480 
481    call timab(545,3,tsec)
482    if ( ((.not.flag_inv_sym) .and. bandpp==1 .and. mpi_enreg%paral_spinor==0 .and. my_nspinor==2 ).or. &
483 &   ((.not.flag_inv_sym) .and. bandpp>1) .or.  flag_inv_sym  ) then
484      if (sij_opt==1) then
485        call xmpi_alltoallv(swavef_alltoall1,recvcountsloc,rdisplsloc,swavef,&
486 &       sendcountsloc,sdisplsloc,spaceComm,ier)
487      end if
488      call xmpi_alltoallv(gvnlc_alltoall1,recvcountsloc,rdisplsloc,gvnlc,&
489 &     sendcountsloc,sdisplsloc,spaceComm,ier)
490      call xmpi_alltoallv(gwavef_alltoall1,recvcountsloc,rdisplsloc,gwavef,&
491 &     sendcountsloc,sdisplsloc,spaceComm,ier)
492    else
493      if (sij_opt==1) then
494        call xmpi_alltoallv(swavef_alltoall2,recvcountsloc,rdisplsloc,swavef,&
495 &       sendcountsloc,sdisplsloc,spaceComm,ier)
496      end if
497      call xmpi_alltoallv(gvnlc_alltoall2,recvcountsloc,rdisplsloc,gvnlc,&
498 &     sendcountsloc,sdisplsloc,spaceComm,ier)
499      call xmpi_alltoallv(gwavef_alltoall2,recvcountsloc,rdisplsloc,gwavef,&
500 &     sendcountsloc,sdisplsloc,spaceComm,ier)
501    end if
502 
503    call timab(545,2,tsec)
504  else
505    if(sij_opt == 1) then
506      call DCOPY(2*ndatarecv*my_nspinor*bandpp, swavef_alltoall2, 1, swavef, 1)
507    end if
508    call DCOPY(2*ndatarecv*my_nspinor*bandpp, gvnlc_alltoall2, 1, gvnlc, 1)
509    call DCOPY(2*ndatarecv*my_nspinor*bandpp, gwavef_alltoall2, 1, gwavef, 1)
510  end if
511 
512 !====================================================================
513  if (flag_inv_sym) then
514    gs_hamk%istwf_k = 2
515  end if
516 !====================================================================
517  ABI_DEALLOCATE(sendcountsloc)
518  ABI_DEALLOCATE(sdisplsloc)
519  ABI_DEALLOCATE(recvcountsloc)
520  ABI_DEALLOCATE(rdisplsloc)
521  ABI_DEALLOCATE(cwavef_alltoall2)
522  ABI_DEALLOCATE(gwavef_alltoall2)
523  ABI_DEALLOCATE(gvnlc_alltoall2)
524  ABI_DEALLOCATE(swavef_alltoall2)
525 
526  if ( ((.not.flag_inv_sym) .and. bandpp==1 .and. mpi_enreg%paral_spinor==0 .and. my_nspinor==2 ).or. &
527 & ((.not.flag_inv_sym) .and. bandpp>1) .or.  flag_inv_sym  ) then
528    ABI_DEALLOCATE(cwavef_alltoall1)
529    ABI_DEALLOCATE(gwavef_alltoall1)
530    ABI_DEALLOCATE(gvnlc_alltoall1)
531    ABI_DEALLOCATE(swavef_alltoall1)
532  end if
533 
534  call timab(630,2,tsec)
535 
536 end subroutine prep_getghc

ABINIT/prep_index_wavef_bandpp [ Functions ]

[ Top ] [ Functions ]

NAME

 prep_index_wavef_bandpp

FUNCTION

 this routine sorts the waves functions by bandpp and by processors
 after the alltoall

INPUTS

  nproc_band = number of processors below the band
  bandpp     = number of groups of couple of waves functions
  nspinor    = number of spin
  ndatarecv  = total number of values received by the processor and sended
               by the other processors band
  recvcounts = number of values sended by each processor band and received
               by the processor
  rdispls    = positions of the values received by the processor and
               sended by each processor band

OUTPUT

  index_wavef_band = position of the sorted values

PARENTS

      chebfi,prep_fourwf,prep_getghc,prep_nonlop

CHILDREN

SOURCE

1940 subroutine prep_index_wavef_bandpp(nproc_band,bandpp,&
1941                              nspinor,ndatarecv,&
1942                              recvcounts,rdispls,&
1943                              index_wavef_band)
1944 
1945 
1946 !This section has been created automatically by the script Abilint (TD).
1947 !Do not modify the following lines by hand.
1948 #undef ABI_FUNC
1949 #define ABI_FUNC 'prep_index_wavef_bandpp'
1950 !End of the abilint section
1951 
1952  implicit none
1953 
1954 !Arguments ------------------------------------
1955 !scalars
1956  integer,intent(in) :: bandpp,ndatarecv,nproc_band,nspinor
1957 !arrays
1958  integer,intent(in) :: rdispls(nproc_band),recvcounts(nproc_band)
1959  integer,allocatable,intent(out) :: index_wavef_band(:)
1960 
1961 !Local variables-------------------------------
1962 !scalars
1963  integer :: delta,idebc,idebe,ifinc,ifine,iindex,iproc,kbandpp,nb
1964 
1965 ! *********************************************************************
1966 
1967 !DEBUG
1968 !write(std_out,*)' prep_index_wavef_banpp : enter '
1969 !write(std_out,*) 'ndatarecv = ', ndatarecv
1970 !write(std_out,*) 'rdispls(:) = ', rdispls(:)
1971 !write(std_out,*) 'recvcounts(:) = ', recvcounts(:)
1972 !ENDDEBUG
1973 
1974 
1975 !---------------------------------------------
1976 !Allocation
1977 !---------------------------------------------
1978  ABI_ALLOCATE(index_wavef_band ,(bandpp*nspinor*ndatarecv))
1979  index_wavef_band(:)   =0
1980 
1981 !---------------------------------------------
1982 !Calcul : loops on bandpp and processors band
1983 !---------------------------------------------
1984  nb = sum(recvcounts(1:nproc_band))
1985  do kbandpp=1,bandpp
1986 
1987    do iproc=1,nproc_band
1988 
1989      idebe = (rdispls(iproc) + 1)  + (kbandpp-1) * ndatarecv*nspinor
1990      ifine = idebe + recvcounts(iproc) -1
1991 
1992      if (iproc==1) then
1993        idebc =   (kbandpp-1)* recvcounts(iproc)*nspinor + 1
1994      else
1995        idebc = (bandpp)  * sum(recvcounts(1:iproc-1))*nspinor &
1996        + (kbandpp-1)* recvcounts(iproc)*nspinor &
1997        + 1
1998      end if
1999      ifinc = idebc + recvcounts(iproc) -1
2000      index_wavef_band(idebe:ifine) = (/( iindex,iindex=idebc,ifinc)/)
2001      delta=ifine-idebe
2002      if (nspinor==2) then
2003        index_wavef_band(idebe+nb :idebe+nb +delta)=(/( iindex,iindex=ifinc+1,ifinc+1+delta)/)
2004      end if
2005    end do
2006  end do
2007 
2008 end subroutine prep_index_wavef_bandpp

abinit/prep_nonlop [ Functions ]

[ Top ] [ abinit ] [ Functions ]

NAME

 prep_nonlop

FUNCTION

 this routine prepares the data to the call of nonlop.

INPUTS

  choice: chooses possible output:
    choice=1 => a non-local energy contribution
          =2 => a gradient with respect to atomic position(s)
          =3 => a gradient with respect to strain(s)
          =23=> a gradient with respect to atm. pos. and strain(s)
          =4 => a 2nd derivative with respect to atomic pos.
          =24=> a gradient and 2nd derivative with respect to atomic pos.
          =5 => a gradient with respect to k wavevector
          =6 => 2nd derivatives with respect to strain and atm. pos.
          =7 => no operator, just projections
  blocksize= size of block for FFT
  cpopt=flag defining the status of cwaveprj=<Proj_i|Cnk> scalars (see below, side effects)
  cwavef(2,npw*my_nspinor*blocksize)=planewave coefficients of wavefunction.
  gvnlc=matrix elements <G|Vnonlocal|C>
  hamk <type(gs_hamiltonian_type)>=data defining the Hamiltonian at a given k (NL part involved here)
  idir=direction of the - atom to be moved in the case (choice=2,signs=2),
                        - k point direction in the case (choice=5,signs=2)
                        - strain component (1:6) in the case (choice=2,signs=2) or (choice=6,signs=1)
  lambdablock(blocksize)=factor to be used when computing (Vln-lambda.S) - only for paw_opt=2
  mpi_enreg=informations about mpi parallelization
  nnlout=dimension of enlout (when signs=1):
  ntypat=number of types of atoms in cell
  paw_opt= define the nonlocal operator concerned with
  signs= if 1, get contracted elements (energy, forces, stress, ...)
         if 2, applies the non-local operator to a function in reciprocal space
  tim_nonlop=timing code of the calling routine (can be set to 0 if not attributed)

OUTPUT

  ==== if (signs==1) ====
  enlout_block(nnlout)=
    if paw_opt==0, 1 or 2: contribution of this block of states to the nl part of various properties
    if paw_opt==3:         contribution of this block of states to <c|S|c>  (where S=overlap when PAW)
  ==== if (signs==2) ====
    if paw_opt==0, 1, 2 or 4:
       gvnlc(2,my_nspinor*npw)=result of the aplication of the nl operator
                        or one of its derivative to the input vect.
    if paw_opt==3 or 4:
       gsc(2,my_nspinor*npw*(paw_opt/3))=result of the aplication of (I+S)
                        to the input vect. (where S=overlap when PAW)

SIDE EFFECTS

  ==== ONLY IF useylm=1
  cwaveprj(natom,my_nspinor) <type(pawcprj_type)>=projected input wave function |c> on non-local projector
                                  =<p_lmn|c> and derivatives
                                  Treatment depends on cpopt parameter:
                     if cpopt=-1, <p_lmn|in> (and derivatives)
                                  have to be computed (and not saved)
                     if cpopt= 0, <p_lmn|in> have to be computed and saved
                                  derivatives are eventually computed but not saved
                     if cpopt= 1, <p_lmn|in> and first derivatives have to be computed and saved
                                  other derivatives are eventually computed but not saved
                     if cpopt= 2  <p_lmn|in> are already in memory;
                                  only derivatives are computed here and not saved
 (if useylm=0, should have cpopt=-1)

PARENTS

      energy,forstrnps,m_invovl,m_lobpcgwf,vtowfk

CHILDREN

      dcopy,nonlop,prep_index_wavef_bandpp,prep_sort_wavef_spin,timab
      xmpi_allgather,xmpi_alltoallv

NOTES

  cprj (as well as cg) is distributed over band processors.
  Only the mod((iband-1)/mpi_enreg%bandpp,mpi_enreg%nproc_band) projected WFs are stored on each proc.

SOURCE

615 subroutine prep_nonlop(choice,cpopt,cwaveprj,enlout_block,hamk,idir,lambdablock,&
616 &                      blocksize,mpi_enreg,nnlout,paw_opt,signs,gsc,&
617 &                      tim_nonlop,cwavef,gvnlc,already_transposed)
618 
619 
620 !This section has been created automatically by the script Abilint (TD).
621 !Do not modify the following lines by hand.
622 #undef ABI_FUNC
623 #define ABI_FUNC 'prep_nonlop'
624 !End of the abilint section
625 
626  implicit none
627 
628 !Arguments ------------------------------------
629  integer,intent(in) :: blocksize,choice,cpopt,idir,signs,nnlout,paw_opt
630  logical,optional,intent(in) :: already_transposed
631  real(dp),intent(in) :: lambdablock(blocksize)
632  real(dp),intent(out) :: enlout_block(nnlout*blocksize),gvnlc(:,:),gsc(:,:)
633  real(dp),intent(inout) :: cwavef(:,:)
634  type(gs_hamiltonian_type),intent(in) :: hamk
635  type(mpi_type),intent(inout) :: mpi_enreg
636  type(pawcprj_type),intent(inout) :: cwaveprj(:,:)
637 
638 !Local variables-------------------------------
639 !scalars
640  integer :: bandpp,ier,ikpt_this_proc,my_nspinor,ndatarecv,nproc_band,npw,nspinortot
641  integer :: old_me_g0,spaceComm=0,tim_nonlop
642  logical :: do_transpose
643  character(len=500) :: msg
644 !arrays
645  integer,allocatable :: index_wavef_band(:)
646  integer,  allocatable :: rdisplsloc(:),recvcountsloc(:),sdisplsloc(:),sendcountsloc(:)
647  integer,ABI_CONTIGUOUS  pointer :: rdispls(:),recvcounts(:),sdispls(:),sendcounts(:)
648  real(dp) :: lambda_nonlop(mpi_enreg%bandpp),tsec(2)
649  real(dp), allocatable :: cwavef_alltoall2(:,:),gvnlc_alltoall2(:,:),gsc_alltoall2(:,:)
650  real(dp), allocatable :: cwavef_alltoall1(:,:),gvnlc_alltoall1(:,:),gsc_alltoall1(:,:)
651  real(dp),allocatable :: enlout(:)
652 
653 ! *************************************************************************
654 
655  DBG_ENTER('COLL')
656 
657  call timab(570,1,tsec)
658 
659  do_transpose = .true.
660  if(present(already_transposed)) then
661    if(already_transposed) do_transpose = .false.
662  end if
663 
664  nproc_band = mpi_enreg%nproc_band
665  bandpp     = mpi_enreg%bandpp
666  spaceComm=mpi_enreg%comm_fft
667  if(mpi_enreg%paral_kgb==1) spaceComm=mpi_enreg%comm_band
668  my_nspinor=max(1,hamk%nspinor/mpi_enreg%nproc_spinor)
669  nspinortot=hamk%nspinor
670 
671 !Check sizes
672  npw=hamk%npw_k;if (.not.do_transpose) npw=hamk%npw_fft_k
673  if (size(cwavef)/=2*npw*my_nspinor*blocksize) then
674    msg = 'Incorrect size for cwavef!'
675    MSG_BUG(msg)
676  end if
677  if(choice/=0.and.signs==2) then
678    if (paw_opt/=3) then
679      if (size(gvnlc)/=2*npw*my_nspinor*blocksize) then
680        msg = 'Incorrect size for gvnlc!'
681        MSG_BUG(msg)
682      end if
683    end if
684    if(paw_opt>=3) then
685      if (size(gsc)/=2*npw*my_nspinor*blocksize) then
686        msg = 'Incorrect size for gsc!'
687        MSG_BUG(msg)
688      end if
689    end if
690  end if
691  if(cpopt>=0) then
692    if (size(cwaveprj)/=hamk%natom*my_nspinor*mpi_enreg%bandpp) then
693      msg = 'Incorrect size for cwaveprj!'
694      MSG_BUG(msg)
695    end if
696  end if
697 
698  ABI_ALLOCATE(sendcountsloc,(nproc_band))
699  ABI_ALLOCATE(sdisplsloc   ,(nproc_band))
700  ABI_ALLOCATE(recvcountsloc,(nproc_band))
701  ABI_ALLOCATE(rdisplsloc   ,(nproc_band))
702 
703  ikpt_this_proc=bandfft_kpt_get_ikpt()
704 
705  recvcounts   => bandfft_kpt(ikpt_this_proc)%recvcounts(:)
706  sendcounts   => bandfft_kpt(ikpt_this_proc)%sendcounts(:)
707  rdispls      => bandfft_kpt(ikpt_this_proc)%rdispls   (:)
708  sdispls      => bandfft_kpt(ikpt_this_proc)%sdispls   (:)
709  ndatarecv    =  bandfft_kpt(ikpt_this_proc)%ndatarecv
710 
711  ABI_ALLOCATE(cwavef_alltoall2,(2,ndatarecv*my_nspinor*bandpp))
712  ABI_ALLOCATE(gsc_alltoall2,(2,ndatarecv*my_nspinor*(paw_opt/3)*bandpp))
713  ABI_ALLOCATE(gvnlc_alltoall2,(2,ndatarecv*my_nspinor*bandpp))
714 
715  if(do_transpose .and. (bandpp/=1 .or. (bandpp==1 .and. mpi_enreg%paral_spinor==0.and.nspinortot==2)))then
716    ABI_ALLOCATE(cwavef_alltoall1,(2,ndatarecv*my_nspinor*bandpp))
717    if(signs==2)then
718      if (paw_opt/=3) then
719        ABI_ALLOCATE(gvnlc_alltoall1,(2,ndatarecv*my_nspinor*bandpp))
720      end if
721      if (paw_opt==3.or.paw_opt==4) then
722        ABI_ALLOCATE(gsc_alltoall1,(2,ndatarecv*my_nspinor*bandpp))
723      end if
724    end if
725  end if
726 
727  ABI_ALLOCATE(enlout,(nnlout*bandpp))
728  enlout = zero
729 
730  recvcountsloc(:)=recvcounts(:)*2*my_nspinor*bandpp
731  rdisplsloc(:)=rdispls(:)*2*my_nspinor*bandpp
732  sendcountsloc(:)=sendcounts(:)*2*my_nspinor
733  sdisplsloc(:)=sdispls(:)*2*my_nspinor
734 
735  if(do_transpose) then
736    call timab(581,1,tsec)
737    if(bandpp/=1 .or. (bandpp==1 .and. mpi_enreg%paral_spinor==0.and.nspinortot==2))then
738      call xmpi_alltoallv(cwavef,sendcountsloc,sdisplsloc,cwavef_alltoall1,&
739 &     recvcountsloc,rdisplsloc,spaceComm,ier)
740    else
741      call xmpi_alltoallv(cwavef,sendcountsloc,sdisplsloc,cwavef_alltoall2,&
742 &     recvcountsloc,rdisplsloc,spaceComm,ier)
743    end if
744    call timab(581,2,tsec)
745  else
746    ! Here, we cheat, and use DCOPY to bypass some compiler's overzealous bound-checking
747    ! (ndatarecv*my_nspinor*bandpp might be greater than the declared size of cwavef)
748    call DCOPY(2*ndatarecv*my_nspinor*bandpp,cwavef,1,cwavef_alltoall2,1)
749  end if
750 
751  if(hamk%istwf_k==2) then
752    old_me_g0=mpi_enreg%me_g0
753    if (mpi_enreg%me_fft==0) then
754      mpi_enreg%me_g0=1
755    else
756      mpi_enreg%me_g0=0
757    end if
758  end if
759 
760 !=====================================================================
761  if (bandpp==1) then
762 
763 
764    if (do_transpose .and. mpi_enreg%paral_spinor==0.and.nspinortot==2) then !Sort WF by spin
765      call prep_sort_wavef_spin(nproc_band,my_nspinor,ndatarecv,recvcounts,rdispls,index_wavef_band)
766      cwavef_alltoall2(:, :)=cwavef_alltoall1(:,index_wavef_band)
767    end if
768 
769    if (paw_opt==2) lambda_nonlop(1)=lambdablock(mpi_enreg%me_band+1)
770    call nonlop(choice,cpopt,cwaveprj,enlout,hamk,idir,lambda_nonlop,mpi_enreg,1,nnlout,paw_opt,&
771 &   signs,gsc_alltoall2,tim_nonlop,cwavef_alltoall2,gvnlc_alltoall2)
772 
773    if (do_transpose .and. mpi_enreg%paral_spinor == 0 .and. nspinortot==2.and.signs==2) then
774      if (paw_opt/=3) gvnlc_alltoall1(:,index_wavef_band)=gvnlc_alltoall2(:,:)
775      if (paw_opt==3.or.paw_opt==4) gsc_alltoall1(:,index_wavef_band)=gsc_alltoall2(:,:)
776    end if
777 
778  else   ! bandpp/=1
779 
780 !  -------------------------------------------------------------
781 !  Computation of the index used to sort the waves functions below bandpp
782 !  -------------------------------------------------------------
783    if(do_transpose) then
784      call prep_index_wavef_bandpp(nproc_band,bandpp,&
785 &     my_nspinor,ndatarecv,recvcounts,rdispls,index_wavef_band)
786 
787 !  -------------------------------------------------------
788 !  Sorting of the waves functions below bandpp
789 !  -------------------------------------------------------
790      cwavef_alltoall2(:,:) = cwavef_alltoall1(:,index_wavef_band)
791    end if
792 
793 !  -------------------------------------------------------
794 !  Call nonlop
795 !  -------------------------------------------------------
796    if(paw_opt == 2) lambda_nonlop(1:bandpp) = lambdablock((mpi_enreg%me_band*bandpp)+1:((mpi_enreg%me_band+1)*bandpp))
797    call nonlop(choice,cpopt,cwaveprj,enlout,hamk,idir,lambda_nonlop,mpi_enreg,bandpp,nnlout,paw_opt,&
798 &   signs,gsc_alltoall2,tim_nonlop,cwavef_alltoall2,gvnlc_alltoall2)
799 
800 !  -----------------------------------------------------
801 !  Sorting of waves functions below the processors
802 !  -----------------------------------------------------
803    if(do_transpose.and.signs==2) then
804      if (paw_opt/=3) gvnlc_alltoall1(:,index_wavef_band)=gvnlc_alltoall2(:,:)
805      if (paw_opt==3.or.paw_opt==4) gsc_alltoall1(:,index_wavef_band)=gsc_alltoall2(:,:)
806    end if
807 
808  end if
809 
810 !=====================================================================
811 !  -------------------------------------------------------
812 !  Deallocation
813 !  -------------------------------------------------------
814  if (allocated(index_wavef_band)) then
815    ABI_DEALLOCATE(index_wavef_band)
816  end if
817 
818 !Transpose the gsc_alltoall or gvlnc_alltoall tabs
819 !according to the paw_opt and signs values
820  if(do_transpose) then
821    if (signs==2) then
822      call timab(581,1,tsec)
823      if(bandpp/=1 .or. (bandpp==1 .and. mpi_enreg%paral_spinor==0.and.nspinortot==2))then
824        if (paw_opt/=3) then
825          call xmpi_alltoallv(gvnlc_alltoall1,recvcountsloc,rdisplsloc,gvnlc,&
826 &         sendcountsloc,sdisplsloc,spaceComm,ier)
827        end if
828        if (paw_opt==3.or.paw_opt==4) then
829          call xmpi_alltoallv(gsc_alltoall1,recvcountsloc,rdisplsloc,gsc,&
830 &         sendcountsloc,sdisplsloc,spaceComm,ier)
831        end if
832      else
833        if (paw_opt/=3) then
834          call xmpi_alltoallv(gvnlc_alltoall2,recvcountsloc,rdisplsloc,gvnlc,&
835 &         sendcountsloc,sdisplsloc,spaceComm,ier)
836        end if
837        if (paw_opt==3.or.paw_opt==4) then
838          call xmpi_alltoallv(gsc_alltoall2,recvcountsloc,rdisplsloc,gsc,&
839 &         sendcountsloc,sdisplsloc,spaceComm,ier)
840        end if
841      end if
842      call timab(581,2,tsec)
843    end if
844  else
845    ! TODO check other usages, maybe
846    call DCOPY(2*ndatarecv*my_nspinor*bandpp, gsc_alltoall2, 1, gsc, 1)
847  end if
848  if (hamk%istwf_k==2) mpi_enreg%me_g0=old_me_g0
849 
850  if (nnlout>0) then
851    call xmpi_allgather(enlout,nnlout*bandpp,enlout_block,spaceComm,ier)
852  end if
853  ABI_DEALLOCATE(enlout)
854  ABI_DEALLOCATE(sendcountsloc)
855  ABI_DEALLOCATE(sdisplsloc)
856  ABI_DEALLOCATE(recvcountsloc)
857  ABI_DEALLOCATE(rdisplsloc)
858  ABI_DEALLOCATE(cwavef_alltoall2)
859  ABI_DEALLOCATE(gvnlc_alltoall2)
860  ABI_DEALLOCATE(gsc_alltoall2)
861  if(do_transpose .and. (bandpp/=1 .or. (bandpp==1 .and. mpi_enreg%paral_spinor==0.and.nspinortot==2)))then
862    ABI_DEALLOCATE(cwavef_alltoall1)
863    if(signs==2)then
864      if (paw_opt/=3) then
865        ABI_DEALLOCATE(gvnlc_alltoall1)
866      end if
867      if (paw_opt==3.or.paw_opt==4) then
868        ABI_DEALLOCATE(gsc_alltoall1)
869      end if
870    end if
871  end if
872 
873  call timab(570,2,tsec)
874 
875  DBG_EXIT('COLL')
876 
877 end subroutine prep_nonlop

ABINIT/prep_sort_wavef_spin [ Functions ]

[ Top ] [ Functions ]

NAME

 prep_sort_wavef_spin

FUNCTION

 Compute index used to sort a spinorial wave-function by spin
 Sort to have all nspinor=1 fisrt, then all nspinor=2

INPUTS

  nproc_band=size of "band" communicator
  nspinor=number of spinorial components of the wavefunction
  ndatarecv=total number of values on all processors
  recvcounts(nproc_band)= number of received values by the processor
  rdispls(nproc_band)= offsets of the received values by the processor

OUTPUT

  index_wavef(:)=array containing the sorted indexes (pointer, allocated in this routine)

PARENTS

      prep_getghc,prep_nonlop

CHILDREN

SOURCE

2036 subroutine prep_sort_wavef_spin(nproc_band,nspinor,ndatarecv,recvcounts,rdispls,index_wavef)
2037 
2038 
2039 !This section has been created automatically by the script Abilint (TD).
2040 !Do not modify the following lines by hand.
2041 #undef ABI_FUNC
2042 #define ABI_FUNC 'prep_sort_wavef_spin'
2043 !End of the abilint section
2044 
2045  implicit none
2046 
2047 !Arguments ------------------------------------
2048 !scalars
2049  integer,intent(in) :: ndatarecv,nproc_band,nspinor
2050 !arrays
2051  integer,intent(in) :: rdispls(nproc_band),recvcounts(nproc_band)
2052  integer,allocatable,intent(out) :: index_wavef(:)
2053 
2054 !Local variables-------------------------------
2055 !scalars
2056  integer :: isft,isft1,iproc,iindex
2057 !arrays
2058  integer,allocatable :: recvcountsloc(:),rdisplsloc(:)
2059 
2060 ! *********************************************************************
2061 
2062  ABI_ALLOCATE(index_wavef,(ndatarecv*nspinor))
2063 
2064  ABI_ALLOCATE(recvcountsloc,(nproc_band))
2065  ABI_ALLOCATE(rdisplsloc,(nproc_band))
2066  recvcountsloc(:)=recvcounts(:)*2*nspinor
2067  rdisplsloc(:)=rdispls(:)*2*nspinor
2068 
2069 !---------------------------------------------
2070 !Loops on bandpp and processors band
2071 !---------------------------------------------
2072  isft=0
2073  do iproc=1,nproc_band
2074 
2075 !  ===== Spin up
2076    if (iproc==1) then
2077      isft= 0
2078    else
2079      isft= sum(recvcounts(1: (iproc-1)))
2080    end if
2081    isft1 = 0.5*rdisplsloc(iproc)
2082 
2083    index_wavef(1+isft:isft+recvcounts(iproc))= &
2084 &   (/(iindex,iindex=isft1+1,isft1+recvcounts(iproc))/)
2085 
2086 !  =====Spin down
2087    if (iproc==1)then
2088      isft=sum(recvcounts(1:nproc_band))
2089    else
2090      isft=sum(recvcounts(1:nproc_band)) &
2091 &     +sum(recvcounts(1:iproc-1))
2092    end if
2093    isft1 = 0.5 * rdisplsloc(iproc) + recvcounts(iproc)
2094 
2095    index_wavef(1+isft:isft+recvcounts(iproc))= &
2096 &   (/(iindex,iindex=isft1+1,isft1+ recvcounts(iproc))/)
2097 
2098  end do
2099 
2100  ABI_DEALLOCATE(recvcountsloc)
2101  ABI_DEALLOCATE(rdisplsloc)
2102 
2103 end subroutine prep_sort_wavef_spin

ABINIT/prep_wavef_sym_do [ Functions ]

[ Top ] [ Functions ]

NAME

 prep_wavef_sym_do

FUNCTION

 this routine associates waves functions by two as following
      E(G)  = C(G) + D(G)
      E(-G) = C*(G) + iD*(G)
 the values are distributed on the processors in function of
 the value of mpi_enreg%distribfft%tab_fftwf2_distrib( (-kg_k_gather(2,i) )

INPUTS

  mpi_enreg          = informations about mpi parallelization
  bandpp             = number of couple of waves functions
  nspinor            = number of spin
  ndatarecv          = number of values received by the processor and sended
                       by the other processors band
  ndatarecv_tot      = total number of received values
                       (ndatarecv   + number of received opposited planewave coordinates)
  ndatasend_sym      = number of sended values to the processors fft to create opposited
                       planewave coordinates
  tab_proc           = positions of opposited planewave coordinates in the list of the
                       processors fft
  cwavef_alltoall    = planewave coefficients of wavefunction
                      ( initial of the processor + sended by other processors band)
  sendcounts_sym     = number of sended values by the processor to each processor fft
  sdispls_sym        = postions of the sended values by the processor to each processor fft

  recvcounts_sym     = number of the received values by the processor from each processor fft
  rdispls_sym        = postions of the received values by the processor from each processor fft

OUTPUT

  ewavef_alltoall_sym = planewave coefficients of wavefunction
                        initial of the processor +
                        sended by other processors band +
                        sended by other processors fft  +
                        and compisited if bandpp >1
  index_wavef_send    = index to send the values in blocks to the other processor fft

SIDE EFFECTS

PARENTS

      prep_fourwf,prep_getghc

CHILDREN

      xmpi_alltoallv

SOURCE

1431 subroutine prep_wavef_sym_do(mpi_enreg,bandpp,nspinor,&
1432 &     ndatarecv,&
1433 &     ndatarecv_tot,ndatasend_sym,tab_proc,&
1434 &     cwavef_alltoall,&
1435 &     sendcounts_sym,sdispls_sym,&
1436 &     recvcounts_sym,rdispls_sym,&
1437 &     ewavef_alltoall_sym,&
1438 &     index_wavef_send)
1439 
1440 
1441 !This section has been created automatically by the script Abilint (TD).
1442 !Do not modify the following lines by hand.
1443 #undef ABI_FUNC
1444 #define ABI_FUNC 'prep_wavef_sym_do'
1445 !End of the abilint section
1446 
1447  implicit none
1448 
1449 !Arguments ------------------------------------
1450 !scalars
1451  integer,intent(in) :: bandpp,ndatarecv,ndatarecv_tot,ndatasend_sym
1452  integer,intent(in) :: nspinor
1453  type(mpi_type),intent(in) :: mpi_enreg
1454 !arrays
1455  integer,allocatable,intent(out) :: index_wavef_send(:)
1456  integer,intent(in) :: rdispls_sym(:),recvcounts_sym(:)
1457  integer,intent(in) :: sdispls_sym(:),sendcounts_sym(:)
1458  integer,intent(in) :: tab_proc(:)
1459  real(dp),intent(inout) :: cwavef_alltoall(2,ndatarecv*nspinor*bandpp)
1460  real(dp),pointer :: ewavef_alltoall_sym(:,:)
1461 
1462 !Local variables-------------------------------
1463 !scalars
1464  integer :: bandpp_sym,ibandpp,idatarecv,ideb_loc,idebc,idebd,idebe
1465  integer :: ier,ifin_loc,ifinc,ifind,ifine,iproc,jbandpp,jsendloc
1466  integer :: kbandpp,newspacecomm,nproc_fft
1467  logical :: flag_compose
1468 !arrays
1469  integer,allocatable :: rdispls_sym_loc(:),recvcounts_sym_loc(:)
1470  integer,allocatable :: sdispls_sym_loc(:),sendcounts_sym_loc(:)
1471  real(dp),allocatable :: ewavef_alltoall_loc(:,:),ewavef_alltoall_send(:,:)
1472 
1473 ! *********************************************************************
1474 
1475 !DEBUG
1476 !write(std_out,*)' prep_wavef_sym_do : enter '
1477 !ENDDEBUG
1478 
1479 !---------------------------------------------
1480 !Initialisation
1481 !---------------------------------------------
1482  nproc_fft    = mpi_enreg%nproc_fft
1483 
1484  newspacecomm = mpi_enreg%comm_fft
1485 
1486  if (modulo(bandpp,2)==0) then
1487    bandpp_sym   = bandpp/2
1488    flag_compose = .TRUE.
1489  else
1490    bandpp_sym   = bandpp
1491    flag_compose = .FALSE.
1492  end if
1493 
1494 !---------------------------------------------
1495 !Allocation
1496 !---------------------------------------------
1497  ABI_ALLOCATE(ewavef_alltoall_sym     ,(2,ndatarecv_tot*bandpp_sym))
1498  ABI_ALLOCATE(ewavef_alltoall_loc     ,(2,ndatarecv    *bandpp_sym))
1499  ABI_ALLOCATE(ewavef_alltoall_send    ,(2,ndatasend_sym*bandpp_sym))
1500  ABI_ALLOCATE(index_wavef_send        ,(  ndatasend_sym*bandpp_sym))
1501 
1502  ABI_ALLOCATE(sendcounts_sym_loc    ,(nproc_fft))
1503  ABI_ALLOCATE(sdispls_sym_loc       ,(nproc_fft))
1504  ABI_ALLOCATE(recvcounts_sym_loc    ,(nproc_fft))
1505  ABI_ALLOCATE(rdispls_sym_loc       ,(nproc_fft))
1506 
1507 
1508 !Initialisation
1509 !--------------
1510  ewavef_alltoall_sym(:,:) =0.
1511  ewavef_alltoall_loc(:,:) =0.
1512 
1513  sendcounts_sym_loc(:) =0
1514  sdispls_sym_loc(:)    =0
1515  recvcounts_sym_loc(:) =0
1516  rdispls_sym_loc(:)    =0
1517 
1518  index_wavef_send(:)   =0
1519 
1520 
1521 !-------------------------------------------------
1522 !We are bandpp blocks which we want to :
1523 !associate by two      (band_sym==bandpp/2)
1524 !or not associate by two  (band_sym==bandpp)
1525 !
1526 !So We'll have got bandpp_sym blocks
1527 !So we loop on the bandpp_sym blocks
1528 !--------------------------------------------------
1529 
1530  do kbandpp=1,bandpp_sym
1531 
1532 !  position of the two blocks
1533 !  --------------------------
1534    ibandpp = (kbandpp-1) * 2
1535    jbandpp =  ibandpp    + 1
1536 
1537    idebe = (kbandpp-1) * ndatarecv_tot + 1
1538    ifine = idebe       + ndatarecv     - 1
1539 
1540    idebc = ibandpp * ndatarecv     + 1
1541    ifinc = idebc   + ndatarecv     - 1
1542 
1543    idebd = jbandpp * ndatarecv     + 1
1544    ifind = idebd   + ndatarecv     - 1
1545 
1546    ideb_loc = (kbandpp-1) * ndatarecv  + 1
1547    ifin_loc = ideb_loc    + ndatarecv  - 1
1548 
1549 
1550    if (flag_compose) then
1551 
1552 
1553 !    calcul ewavef(G)
1554 !    ----------------
1555      ewavef_alltoall_sym(1,idebe:ifine) =    &
1556 &     cwavef_alltoall(1,idebc:ifinc) &
1557 &     - cwavef_alltoall(2,idebd:ifind)
1558 
1559      ewavef_alltoall_sym(2,idebe:ifine) =    &
1560 &     cwavef_alltoall(2,idebc:ifinc) &
1561 &     + cwavef_alltoall(1,idebd:ifind)
1562 
1563 !    calcul ewavef_loc(-G)
1564 !    ---------------------
1565      ewavef_alltoall_loc(1,ideb_loc:ifin_loc) =  &
1566 &     cwavef_alltoall(1,idebc:ifinc) &
1567 &     + cwavef_alltoall(2,idebd:ifind)
1568 
1569      ewavef_alltoall_loc(2,ideb_loc:ifin_loc) =  &
1570 &     - cwavef_alltoall(2,idebc:ifinc) &
1571 &     + cwavef_alltoall(1,idebd:ifind)
1572    else
1573 
1574 !    calcul ewavef(G)
1575 !    ----------------
1576      ewavef_alltoall_sym(1,idebe:ifine)   = cwavef_alltoall(1,idebc:ifinc)
1577      ewavef_alltoall_sym(2,idebe:ifine)   = cwavef_alltoall(2,idebc:ifinc)
1578 
1579 !    calcul ewavef_loc(-G)
1580 !    ---------------------
1581      ewavef_alltoall_loc(1,ideb_loc:ifin_loc) =   cwavef_alltoall(1,idebc:ifinc)
1582      ewavef_alltoall_loc(2,ideb_loc:ifin_loc) = - cwavef_alltoall(2,idebc:ifinc)
1583 
1584    end if
1585 
1586  end do
1587 
1588 
1589 
1590 !------------------------------------------------------------------------
1591 !Creation of datas blocks for each processor fft from ewavef_alltoall_loc
1592 !to send datas by blocks with a alltoall...
1593 !------------------------------------------------------------------------
1594 
1595 !Position of the blocks
1596 !----------------------
1597  jsendloc=0
1598  do ibandpp=1,bandpp_sym
1599    do iproc=1,nproc_fft
1600      do idatarecv=1,ndatarecv
1601        if (tab_proc(idatarecv)==(iproc-1)) then
1602          jsendloc=jsendloc+1
1603          index_wavef_send(jsendloc)  = idatarecv + ndatarecv * (ibandpp-1)
1604        end if
1605      end do
1606    end do
1607  end do
1608 
1609 !Classment
1610 !----------
1611  ewavef_alltoall_send(:,:)=ewavef_alltoall_loc(:,index_wavef_send)
1612 
1613 
1614 !-------------------------------------------------
1615 !Calcul of the number of received and sended datas
1616 !-------------------------------------------------
1617  sendcounts_sym_loc = sendcounts_sym*2
1618  recvcounts_sym_loc = recvcounts_sym*2
1619 
1620 !------------------------------------------
1621 !Exchange of the datas ewavef_allto_all_loc
1622 !------------------------------------------
1623  do ibandpp=1,bandpp_sym
1624 
1625 !  ------------------------------------------------
1626 !  Deplacment of the sended datas because of bandpp
1627 !  ------------------------------------------------
1628    sdispls_sym_loc(:) = sdispls_sym(:) + ndatasend_sym * (ibandpp-1)
1629    sdispls_sym_loc    = sdispls_sym_loc   *2
1630 
1631 !  --------------------------------------------------
1632 !  Deplacment of the received datas because of bandpp
1633 !  --------------------------------------------------
1634    rdispls_sym_loc(:) = rdispls_sym(:) + ndatarecv_tot * (ibandpp-1)
1635    rdispls_sym_loc    = rdispls_sym_loc   *2
1636 
1637 
1638    call xmpi_alltoallv(&
1639 &   ewavef_alltoall_send(:,:) ,sendcounts_sym_loc,sdispls_sym_loc,&
1640 &   ewavef_alltoall_sym(:,:)  ,recvcounts_sym_loc,rdispls_sym_loc,&
1641 &   newspacecomm,ier)
1642 
1643  end do
1644 
1645 !-----------------------
1646 !Desallocation
1647 !-----------------------
1648 
1649  ABI_DEALLOCATE(sendcounts_sym_loc)
1650  ABI_DEALLOCATE(recvcounts_sym_loc)
1651  ABI_DEALLOCATE(sdispls_sym_loc)
1652  ABI_DEALLOCATE(rdispls_sym_loc)
1653 
1654  ABI_DEALLOCATE(ewavef_alltoall_loc)
1655  ABI_DEALLOCATE(ewavef_alltoall_send)
1656 
1657 end subroutine prep_wavef_sym_do

ABINIT/prep_wavef_sym_undo [ Functions ]

[ Top ] [ Functions ]

NAME

 prep_wavef_sym_undo

FUNCTION

 this routine dissociates each wave function in two waves functions as following
      C(G) =   ( E*(-G) + E(G))/2
      D(G) = i*( E*(-G) - E(G))/2
 the values are redistributed on the processors in function of
 the value of mpi_enreg%distribfft%tab_fftwf2_distrib( (-kg_k_gather(2,i) )

INPUTS

  mpi_enreg          = informations about mpi parallelization
  bandpp             = number of groups of couple of waves functions
  nspinor            = number of spin
  ndatarecv          = number of values received by the processor and sended
                       by the other processors band
  ndatarecv_tot      = total number of received values
                       (ndatarecv   + number of received opposited planewave coordinates)
  ndatasend_sym      = number of sended values to the processors fft to create opposited
                       planewave coordinates
  idatarecv0         = position of the planewave coordinates (0,0,0)
  sendcounts_sym     = number of sended values by the processor to each processor fft
  sdispls_sym        = postions of the sended values by the processor to each processor fft

  recvcounts_sym     = number of the received values by the processor to each processor fft

  gwavef_alltoall_sym = planewave coefficients of wavefunction
                        initial of the processor +
                        sended by other processors band +
                        sended by other processors fft  +
                        and composited if bandpp >1
  index_wavef_send    = index to send the values by block to the other processor fft

OUTPUT

  gwavef_alltoall     = planewave coefficients of wavefunction
                        ( for of the processor + to send to other processors band)

PARENTS

      prep_fourwf,prep_getghc

CHILDREN

      xmpi_alltoallv

SOURCE

1707 subroutine prep_wavef_sym_undo(mpi_enreg,bandpp,nspinor,&
1708 &     ndatarecv,&
1709 &     ndatarecv_tot,ndatasend_sym,idatarecv0,&
1710 &     gwavef_alltoall,&
1711 &     sendcounts_sym,sdispls_sym,&
1712 &     recvcounts_sym,rdispls_sym,&
1713 &     gwavef_alltoall_sym,&
1714 &     index_wavef_send)
1715 
1716 
1717 !This section has been created automatically by the script Abilint (TD).
1718 !Do not modify the following lines by hand.
1719 #undef ABI_FUNC
1720 #define ABI_FUNC 'prep_wavef_sym_undo'
1721 !End of the abilint section
1722 
1723  implicit none
1724 
1725 !Arguments ------------------------------------
1726 !scalars
1727  integer,intent(in) :: bandpp,idatarecv0,ndatarecv,ndatarecv_tot,ndatasend_sym
1728  integer,intent(in) :: nspinor
1729  type(mpi_type),intent(in) :: mpi_enreg
1730 !arrays
1731  integer,intent(in) :: index_wavef_send(:),rdispls_sym(:),recvcounts_sym(:)
1732  integer,intent(in) :: sdispls_sym(:),sendcounts_sym(:)
1733  real(dp),intent(inout) :: gwavef_alltoall(2,ndatarecv*nspinor*bandpp)
1734  real(dp),intent(inout) :: gwavef_alltoall_sym(:,:)
1735 
1736 !Local variables-------------------------------
1737 !scalars
1738  integer :: bandpp_sym,ibandpp,ideb_loc,idebc,idebd
1739  integer :: idebe,ier,ifin_loc,ifinc,ifind,ifine
1740  integer :: jbandpp,kbandpp,newspacecomm,nproc_fft
1741  logical :: flag_compose
1742 !arrays
1743  integer,allocatable :: rdispls_sym_loc(:),recvcounts_sym_loc(:)
1744  integer,allocatable :: sdispls_sym_loc(:),sendcounts_sym_loc(:)
1745  real(dp),allocatable :: gwavef_alltoall_loc(:,:),gwavef_alltoall_rcv(:,:)
1746 
1747 ! *********************************************************************
1748 
1749 !DEBUG
1750 !write(std_out,*)' prep_wavef_sym_undo : enter '
1751 !ENDDEBUG
1752 
1753 
1754 !---------------------------------------------
1755 !Initialisation
1756 !---------------------------------------------
1757  nproc_fft    = mpi_enreg%nproc_fft
1758 
1759  newspacecomm = mpi_enreg%comm_fft
1760 
1761  if (modulo(bandpp,2)==0) then
1762    bandpp_sym   = bandpp/2
1763    flag_compose = .TRUE.
1764  else
1765    bandpp_sym   = bandpp
1766    flag_compose = .FALSE.
1767  end if
1768 
1769 !---------------------------------------------
1770 !Allocation
1771 !---------------------------------------------
1772  ABI_ALLOCATE(gwavef_alltoall_loc     ,(2,ndatarecv     *bandpp_sym))
1773  ABI_ALLOCATE(gwavef_alltoall_rcv     ,(2,ndatasend_sym *bandpp_sym))
1774 
1775  ABI_ALLOCATE(sendcounts_sym_loc    ,(nproc_fft))
1776  ABI_ALLOCATE(sdispls_sym_loc       ,(nproc_fft))
1777  ABI_ALLOCATE(recvcounts_sym_loc    ,(nproc_fft))
1778  ABI_ALLOCATE(rdispls_sym_loc       ,(nproc_fft))
1779 
1780 
1781 !---------------------------------------------
1782 !Initialisation
1783 !---------------------------------------------
1784  gwavef_alltoall_loc(:,:) =0.
1785 
1786  sendcounts_sym_loc(:) =0
1787  sdispls_sym_loc(:)    =0
1788  recvcounts_sym_loc(:) =0
1789  rdispls_sym_loc(:)    =0
1790 
1791 
1792 !-------------------------------------------------
1793 !Calcul of number of the sended and received datas
1794 !-------------------------------------------------
1795  sendcounts_sym_loc = sendcounts_sym*2
1796  recvcounts_sym_loc = recvcounts_sym*2
1797 
1798 !----------------------------------------------------
1799 !Sending of the values
1800 !----------------------------------------------------
1801  do ibandpp = 1,bandpp_sym
1802 
1803 !  -------------------------------------------------
1804 !  Deplacment of the sended values because of bandpp
1805 !  -------------------------------------------------
1806    sdispls_sym_loc(:) = sdispls_sym(:) + ndatasend_sym * (ibandpp-1)
1807    sdispls_sym_loc    = sdispls_sym_loc   *2
1808 
1809 !  ---------------------------------------------------
1810 !  Deplacment of the received values because of bandpp
1811 !  ---------------------------------------------------
1812    rdispls_sym_loc(:) = rdispls_sym(:) + ndatarecv_tot * (ibandpp-1)
1813    rdispls_sym_loc    = rdispls_sym_loc   *2
1814 
1815 
1816    call xmpi_alltoallv(&
1817 &   gwavef_alltoall_sym(:,:) ,recvcounts_sym_loc,rdispls_sym_loc,&
1818 &   gwavef_alltoall_rcv(:,:) ,sendcounts_sym_loc,sdispls_sym_loc,&
1819 &   newspacecomm,ier)
1820 
1821  end do
1822 
1823 
1824 !----------------------
1825 !Dispatching the blocks
1826 !----------------------
1827  gwavef_alltoall_loc(:,index_wavef_send(:)) = gwavef_alltoall_rcv(:,:)
1828 
1829 !----------------------
1830 !Case  -kg = [0 0 0]
1831 !----------------------
1832  if (idatarecv0/=-1) then
1833    do kbandpp=1,bandpp_sym
1834      gwavef_alltoall_loc(:,(kbandpp-1)*ndatarecv     + idatarecv0)= &
1835      gwavef_alltoall_sym(:,(kbandpp-1)*ndatarecv_tot + idatarecv0)
1836    end do
1837  end if
1838 
1839 !---------------------------------------------------
1840 !Build of hwavef_alltoall
1841 !
1842 !We have got :
1843 !bandpp_sym blocks to dissociate
1844 !or   bandpp_sym blokcs to not dissociate
1845 !--------------------------------------------------
1846  do kbandpp=1,bandpp_sym
1847 
1848 !  position of the 2 blocks
1849 !  ----------------------------------
1850    ibandpp = (kbandpp-1) * 2
1851    jbandpp =  ibandpp    + 1
1852 
1853    idebe = (kbandpp-1) * ndatarecv_tot + 1
1854    ifine = idebe       + ndatarecv     - 1
1855 
1856    idebc = ibandpp * ndatarecv     + 1
1857    ifinc = idebc   + ndatarecv     - 1
1858 
1859    idebd = jbandpp * ndatarecv     + 1
1860    ifind = idebd   + ndatarecv     - 1
1861 
1862    ideb_loc = (kbandpp-1) * ndatarecv  + 1
1863    ifin_loc = ideb_loc    + ndatarecv  - 1
1864 
1865 
1866    if (flag_compose) then
1867 
1868 !    calcul cwavef(G)
1869 !    ----------------
1870      gwavef_alltoall(1,idebc:ifinc) =   gwavef_alltoall_sym(1,idebe:ifine)  &
1871 &     + gwavef_alltoall_loc(1,ideb_loc:ifin_loc)
1872      gwavef_alltoall(2,idebc:ifinc) =   gwavef_alltoall_sym(2,idebe:ifine)  &
1873 &     - gwavef_alltoall_loc(2,ideb_loc:ifin_loc)
1874 
1875 !    calcul dwavef(G)
1876 !    ------------------
1877      gwavef_alltoall(1,idebd:ifind) =   gwavef_alltoall_sym(2,idebe:ifine) &
1878 &     + gwavef_alltoall_loc(2,ideb_loc:ifin_loc)
1879      gwavef_alltoall(2,idebd:ifind) = - gwavef_alltoall_sym(1,idebe:ifine) &
1880 &     + gwavef_alltoall_loc(1,ideb_loc:ifin_loc)
1881    else
1882 
1883 !    calcul cwavef(G)
1884 !    ----------------
1885      gwavef_alltoall(1,idebc:ifinc) =   gwavef_alltoall_sym(1,idebe:ifine)  &
1886 &     + gwavef_alltoall_loc(1,ideb_loc:ifin_loc)
1887      gwavef_alltoall(2,idebc:ifinc) =   gwavef_alltoall_sym(2,idebe:ifine)  &
1888 &     - gwavef_alltoall_loc(2,ideb_loc:ifin_loc)
1889    end if
1890 
1891  end do
1892 
1893 !We divise by two
1894  gwavef_alltoall(:,:)    = gwavef_alltoall(:,:)/2
1895 
1896 !-----------------------
1897 !Desallocation
1898 !-----------------------
1899 
1900  ABI_DEALLOCATE(sendcounts_sym_loc)
1901  ABI_DEALLOCATE(recvcounts_sym_loc)
1902  ABI_DEALLOCATE(sdispls_sym_loc)
1903  ABI_DEALLOCATE(rdispls_sym_loc)
1904 
1905  ABI_DEALLOCATE(gwavef_alltoall_loc)
1906  ABI_DEALLOCATE(gwavef_alltoall_rcv)
1907 
1908 end subroutine prep_wavef_sym_undo