TABLE OF CONTENTS


ABINIT/m_wfk_analyze [ Modules ]

[ Top ] [ Modules ]

NAME

  m_wfk_analyze

FUNCTION

  Post-processing tools for WFK file

COPYRIGHT

  Copyright (C) 2008-2018 ABINIT group (MG)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_wfk_analyze
28 
29  use defs_basis
30  use m_errors
31  use m_abicore
32 
33  implicit none
34 
35  private

ABINIT/wfk_analyze [ Functions ]

[ Top ] [ Functions ]

NAME

  wfk_analyze

FUNCTION

 Main routine implementing postprocessing tools for the WFK file.
 Main differences wrt cut3d:

   - MPI support.
   - No interactive prompt.
   - Run the analysis once, store all the important results in netcdf files
     and use python tools to analyze data

INPUTS

 acell(3)=Length scales of primitive translations (bohr)
 codvsn=Code version
 dtfil<datafiles_type>=Variables related to files.
 dtset<dataset_type>=All input variables for this dataset.
 pawang<pawang_type)>=PAW angular mesh and related data.
 pawrad(ntypat*usepaw)<pawrad_type>=Paw radial mesh and related data.
 pawtab(ntypat*usepaw)<pawtab_type>=Paw tabulated starting data.
 psps<pseudopotential_type>=Variables related to pseudopotentials.
   Before entering the first time in the routine, a significant part of psps has been initialized :
   the integers dimekb,lmnmax,lnmax,mpssang,mpssoang,mpsso,mgrid,ntypat,n1xccc,usepaw,useylm,
   and the arrays dimensioned to npsp. All the remaining components of psps are to be initialized in
   the call to pspini. The next time the code enters bethe_salpeter, psps might be identical to the
   one of the previous dtset, in which case, no reinitialisation is scheduled in pspini.F90.
 rprim(3,3)=Dimensionless real space primitive translations.
 xred(3,natom)=Reduced atomic coordinates.

PARENTS

      driver

NOTES

 ON THE USE OF FFT GRIDS:
 =================
 In case of PAW:
 ---------------
    Two FFT grids are used:
    - A "coarse" FFT grid (defined by ecut) for the application of the Hamiltonian on the plane waves basis.
      It is defined by nfft, ngfft, mgfft, ...
      Hamiltonian, wave-functions, density related to WFs (rhor here), ... are expressed on this grid.
    - A "fine" FFT grid (defined) by ecutdg) for the computation of the density inside PAW spheres.
      It is defined by nfftf, ngfftf, mgfftf, ... Total density, potentials, ... are expressed on this grid.
 In case of norm-conserving:
 ---------------------------
    - Only the usual FFT grid (defined by ecut) is used. It is defined by nfft, ngfft, mgfft, ...
      For compatibility reasons, (nfftf,ngfftf,mgfftf) are set equal to (nfft,ngfft,mgfft) in that case.

CHILDREN

      wfd_init,wfd_read_wfk,wfd_test_ortho

SOURCE

 99 subroutine wfk_analyze(acell,codvsn,dtfil,dtset,pawang,pawrad,pawtab,psps,rprim,xred)
100 
101  use defs_basis
102  use defs_datatypes
103  use defs_abitypes
104  use m_abicore
105  use m_xmpi
106  use m_errors
107  use m_hdr
108  use m_crystal
109  use m_crystal_io
110  use m_ebands
111  use m_tetrahedron
112  use m_nctk
113  use m_wfk
114  use m_wfd
115 
116  use m_time,            only : timab
117  use m_fstrings,        only : strcat, sjoin, itoa, ftoa
118  use m_fftcore,         only : print_ngfft
119  use m_kpts,            only : tetra_from_kptrlatt
120  use m_bz_mesh,         only : kpath_t, kpath_new, kpath_free
121  use m_mpinfo,          only : destroy_mpi_enreg, initmpi_seq
122  use m_esymm,           only : esymm_t, esymm_free, esymm_failed
123  use m_ddk,             only : eph_ddk
124  use m_pawang,          only : pawang_type
125  use m_pawrad,          only : pawrad_type
126  use m_pawtab,          only : pawtab_type, pawtab_print, pawtab_get_lsize
127  use m_paw_an,          only : paw_an_type, paw_an_init, paw_an_free, paw_an_nullify
128  use m_paw_ij,          only : paw_ij_type, paw_ij_init, paw_ij_free, paw_ij_nullify
129  use m_pawfgrtab,       only : pawfgrtab_type, pawfgrtab_free, pawfgrtab_init, pawfgrtab_print
130  use m_pawrhoij,        only : pawrhoij_type, pawrhoij_alloc, pawrhoij_copy, pawrhoij_free, pawrhoij_get_nspden, symrhoij
131  use m_pawdij,          only : pawdij, symdij
132  use m_pawfgr,          only : pawfgr_type, pawfgr_init, pawfgr_destroy
133  use m_paw_sphharm,     only : setsym_ylm
134  use m_paw_init,        only : pawinit,paw_gencond
135  use m_paw_nhat,        only : nhatgrid
136  use m_paw_tools,       only : chkpawovlp
137  use m_paw_correlations,only : pawpuxinit
138  !use m_pawpwij,        only : pawpwff_t, pawpwff_init, pawpwff_free
139  !use m_paw_dmft,       only : paw_dmft_type
140  use m_paw_pwaves_lmn,  only : paw_pwaves_lmn_t, paw_pwaves_lmn_init, paw_pwaves_lmn_free
141  use m_classify_bands,  only : classify_bands
142  use m_pspini,          only : pspini
143 
144 !This section has been created automatically by the script Abilint (TD).
145 !Do not modify the following lines by hand.
146 #undef ABI_FUNC
147 #define ABI_FUNC 'wfk_analyze'
148 !End of the abilint section
149 
150  implicit none
151 
152 !Arguments ------------------------------------
153 !scalars
154  character(len=6),intent(in) :: codvsn
155  type(datafiles_type),intent(in) :: dtfil
156  type(dataset_type),intent(in) :: dtset
157  type(pawang_type),intent(inout) :: pawang
158  type(pseudopotential_type),intent(inout) :: psps
159 !arrays
160  real(dp),intent(in) :: acell(3),rprim(3,3),xred(3,dtset%natom)
161  type(pawrad_type),intent(inout) :: pawrad(psps%ntypat*psps%usepaw)
162  type(pawtab_type),intent(inout) :: pawtab(psps%ntypat*psps%usepaw)
163 
164 !Local variables ------------------------------
165 !scalars
166  integer,parameter :: master=0,brav1=1,timrev2=2,dummy_npw=1
167  integer :: comm,nprocs,my_rank,mgfftf,nfftf !,nfftf_tot
168  integer :: optcut,optgr0,optgr1,optgr2,optrad,psp_gencond !,ii
169  !integer :: option,option_test,option_dij,optrhoij
170  integer :: band,ik_ibz,spin,first_band,last_band
171  integer :: ierr,usexcnhat !,edos_intmeth
172  integer :: cplex_dij,cplex,ndij,nspden_rhoij,gnt_option
173  real(dp),parameter :: spinmagntarget=-99.99_dp
174  real(dp) :: ecore,ecut_eff,ecutdg_eff,gsqcutc_eff,gsqcutf_eff,gsqcut_shp
175  !real(dp) :: edos_step,edos_broad
176  !real(dp) :: cpu,wall,gflops
177  !real(dp) :: ex_energy,gsqcutc_eff,gsqcutf_eff,nelect,norm,oldefermi
178  character(len=500) :: msg
179  character(len=fnlen) :: wfk0_path,wfkfull_path
180  logical :: call_pawinit, use_paw_aeur
181  type(hdr_type) :: wfk0_hdr
182  type(crystal_t) :: cryst
183  type(ebands_t) :: ebands
184  type(t_tetrahedron) :: tetra
185  !type(edos_t) :: edos
186  type(pawfgr_type) :: pawfgr
187  !type(paw_dmft_type) :: paw_dmft
188  type(mpi_type) :: mpi_enreg
189  type(wfd_t) :: wfd
190 !arrays
191  integer :: dummy_gvec(3,dummy_npw),ngfftc(18),ngfftf(18)
192  integer,allocatable :: l_size_atm(:)
193  real(dp),parameter :: k0(3)=zero
194  !real(dp) :: nelect_per_spin(dtset%nsppol),n0(dtset%nsppol)
195  real(dp),pointer :: gs_eigen(:,:,:)
196  complex(gwpc),allocatable :: ur_ae(:)
197  logical,allocatable :: keep_ur(:,:,:),bks_mask(:,:,:)
198  real(dp) :: tsec(2)
199  type(Pawrhoij_type),allocatable :: pawrhoij(:)
200  type(pawfgrtab_type),allocatable :: pawfgrtab(:)
201  !type(paw_ij_type),allocatable :: paw_ij(:)
202  !type(paw_an_type),allocatable :: paw_an(:)
203  type(esymm_t),allocatable :: esymm(:,:)
204  type(paw_pwaves_lmn_t),allocatable :: Paw_onsite(:)
205 
206 !************************************************************************
207 
208  DBG_ENTER('COLL')
209 
210  ! abirules!
211  if (.False.) write(std_out,*)acell,codvsn,rprim,xred
212 
213  comm = xmpi_world; nprocs = xmpi_comm_size(comm); my_rank = xmpi_comm_rank(comm)
214 
215  wfk0_path = dtfil%fnamewffk
216  if (my_rank == master) then
217    ! Accept WFK file in Fortran or netcdf format.
218    if (nctk_try_fort_or_ncfile(wfk0_path, msg) /= 0) then
219      MSG_ERROR(sjoin("Cannot find GS WFK file:", ch10, msg))
220    end if
221  end if
222  call xmpi_bcast(wfk0_path,master,comm,ierr)
223  call wrtout(ab_out, sjoin("- Reading GS states from WFK file:", wfk0_path) )
224 
225  !call cwtime(cpu,wall,gflops,"start")
226 
227  ! Costruct crystal and ebands from the GS WFK file.
228  call wfk_read_eigenvalues(wfk0_path,gs_eigen,wfk0_hdr,comm) !,gs_occ)
229  call hdr_vs_dtset(wfk0_hdr,dtset)
230 
231  call crystal_from_hdr(cryst,wfk0_hdr,timrev2)
232  call crystal_print(cryst,header="crystal structure from WFK file")
233 
234  ebands = ebands_from_hdr(wfk0_hdr,maxval(wfk0_hdr%nband),gs_eigen)
235 
236  ! TODO:
237  ! Make sure everything is OK if WFK comes from a NSCF run since occ are set to zero
238  ! fermie is set to 0 if nscf!
239 
240  ! Here we change the GS bands (fermi level, scissors operator ...)
241  ! All the modifications to ebands should be done here.
242 
243  !if (dtset%occopt /= ebands%occopt .or. abs(dtset%tsmear - ebands%tsmear) > tol12) then
244  !  write(msg,"(2a,2(a,i0,a,f14.6,a))")&
245  !    " Changing occupation scheme as input occopt and tsmear differ from those read from WFK file.",ch10,&
246  !    "   From WFK file: occopt = ",ebands%occopt,", tsmear = ",ebands%tsmear,ch10,&
247  !    "   From input:    occopt = ",dtset%occopt,", tsmear = ",dtset%tsmear,ch10
248  !  call wrtout(ab_out,msg)
249  !  call ebands_set_scheme(ebands,dtset%occopt,dtset%tsmear,spinmagntarget,dtset%prtvol)
250  !end if
251 
252  !if (dtset%eph_fermie /= zero) then ! default value of eph_fermie is zero hence no tolerance is used!
253  !  ABI_CHECK(abs(dtset%eph_extrael) <= tol12, "eph_fermie and eph_extrael are mutually exclusive")
254  !  call wrtout(ab_out, sjoin(" Fermi level set by the user at:",ftoa(dtset%eph_fermie)))
255  !  call ebands_set_fermie(ebands, dtset%eph_fermie, msg)
256  !  call wrtout(ab_out,msg)
257 
258  !else if (abs(dtset%eph_extrael) > tol12) then
259  !  NOT_IMPLEMENTED_ERROR()
260  !  ! TODO: Be careful with the trick used in elphon for passing the concentration
261  !  !call ebands_set_nelect(ebands, dtset%eph_extrael, spinmagntarget, msg)
262  !  call wrtout(ab_out,msg)
263  !end if
264 
265  !call ebands_update_occ(ebands, spinmagntarget)
266  call ebands_print(ebands,header="Ground state energies",prtvol=dtset%prtvol)
267  ABI_FREE(gs_eigen)
268 
269  ! Compute electron DOS.
270  ! TODO: Optimize this part. Really slow if tetra and lots of points
271  ! Could just do DOS around efermi
272  !edos_intmeth = 2; if (dtset%prtdos == 1) edos_intmeth = 1
273  !edos_step = dtset%dosdeltae; edos_broad = dtset%tsmear
274  !edos_step = 0.01 * eV_Ha; edos_broad = 0.3 * eV_Ha
275  !call edos_init(edos,ebands,cryst,edos_intmeth,edos_step,edos_broad,comm,ierr)
276  !ABI_CHECK(ierr==0, "Error in edos_init, see message above.")
277 
278  ! Store DOS per spin channels
279  !n0(:) = edos%gef(1:edos%nsppol)
280  !if (my_rank == master) then
281  !  path = strcat(dtfil%filnam_ds(4), "_EDOS")
282  !  call edos_write(edos, path)
283  !  !call edos_print(edos)
284  !  write(ab_out,"(a)")sjoin("- Writing electron DOS to file:", path)
285  !  write(ab_out,'(a,es16.8,a)')' Fermi level: ',edos%mesh(edos%ief)*Ha_eV," [eV]"
286  !  write(ab_out,"(a,es16.8)")" Total electron DOS in states/eV : ",edos%gef(0) / Ha_eV
287  !  if (ebands%nsppol == 2) then
288  !    write(ab_out,"(a,es16.8)")"   Spin up:  ",edos%gef(1) / Ha_eV
289  !    write(ab_out,"(a,es16.8)")"   Spin down:",edos%gef(2) / Ha_eV
290  !  end if
291  !end if
292  !call edos_free(edos)
293 
294  ! Output useful info on the electronic bands.
295  ! Fermi Surface
296  !if (dtset%prtfsurf /= 0  .and. my_rank == master) then
297  !  path = strcat(dtfil%filnam_ds(4), "_BXSF")
298  !  if (ebands_write_bxsf(ebands,cryst,path) /= 0) then
299  !    MSG_WARNING("Cannot produce file for Fermi surface, check log file for more info")
300  !  end if
301  !end if
302 
303  ! TODO Recheck getng, should use same trick as that used in screening and sigma.
304  call pawfgr_init(pawfgr,dtset,mgfftf,nfftf,ecut_eff,ecutdg_eff,ngfftc,ngfftf,&
305 & gsqcutc_eff=gsqcutc_eff,gsqcutf_eff=gsqcutf_eff,gmet=cryst%gmet,k0=k0)
306 
307  call print_ngfft(ngfftc,header='Coarse FFT mesh used for the wavefunctions')
308  call print_ngfft(ngfftf,header='Dense FFT mesh used for densities and potentials')
309 
310  ! Fake MPI_type for the sequential part.
311  call initmpi_seq(mpi_enreg)
312  call init_distribfft_seq(mpi_enreg%distribfft,'c',ngfftc(2),ngfftc(3),'all')
313  call init_distribfft_seq(mpi_enreg%distribfft,'f',ngfftf(2),ngfftf(3),'all')
314 
315  ! ===========================================
316  ! === Open and read pseudopotential files ===
317  ! ===========================================
318  call pspini(dtset,dtfil,ecore,psp_gencond,gsqcutc_eff,gsqcutf_eff,&
319 & pawrad,pawtab,psps,cryst%rprimd,comm_mpi=comm)
320 
321  ! ============================
322  ! ==== PAW initialization ====
323  ! ============================
324  if (dtset%usepaw == 1) then
325    call chkpawovlp(cryst%natom,cryst%ntypat,dtset%pawovlp,pawtab,cryst%rmet,cryst%typat,cryst%xred)
326 
327    cplex_dij=dtset%nspinor; cplex=1; ndij=1
328 
329    ABI_DT_MALLOC(pawrhoij,(cryst%natom))
330    nspden_rhoij=pawrhoij_get_nspden(dtset%nspden,dtset%nspinor,dtset%pawspnorb)
331    call pawrhoij_alloc(pawrhoij,dtset%pawcpxocc,nspden_rhoij,dtset%nspinor,dtset%nsppol,cryst%typat,pawtab=pawtab)
332 
333    ! Initialize values for several basic arrays
334    gnt_option=1;if (dtset%pawxcdev==2.or.(dtset%pawxcdev==1.and.dtset%positron/=0)) gnt_option=2
335 
336    ! Test if we have to call pawinit
337    call paw_gencond(dtset,gnt_option,"test",call_pawinit)
338 
339    if (psp_gencond==1 .or. call_pawinit) then
340      call timab(553,1,tsec)
341      gsqcut_shp = two*abs(dtset%diecut)*dtset%dilatmx**2/pi**2
342      call pawinit(gnt_option,gsqcut_shp,zero,dtset%pawlcutd,dtset%pawlmix,&
343 &     psps%mpsang,dtset%pawnphi,cryst%nsym,dtset%pawntheta,pawang,Pawrad,&
344 &     dtset%pawspnorb,pawtab,dtset%pawxcdev,dtset%xclevel,dtset%usepotzero)
345      call timab(553,2,tsec)
346 
347      ! Update internal values
348      call paw_gencond(dtset,gnt_option,"save",call_pawinit)
349 
350    else
351      if (pawtab(1)%has_kij  ==1) pawtab(1:cryst%ntypat)%has_kij  =2
352      if (pawtab(1)%has_nabla==1) pawtab(1:cryst%ntypat)%has_nabla=2
353    end if
354 
355    psps%n1xccc=MAXVAL(pawtab(1:cryst%ntypat)%usetcore)
356 
357    ! Initialize optional flags in pawtab to zero
358    ! (Cannot be done in Pawinit since the routine is called only if some pars. are changed)
359    pawtab(:)%has_nabla = 0
360    pawtab(:)%usepawu   = 0
361    pawtab(:)%useexexch = 0
362    pawtab(:)%exchmix   =zero
363 
364    call setsym_ylm(cryst%gprimd,pawang%l_max-1,cryst%nsym,dtset%pawprtvol,cryst%rprimd,cryst%symrec,pawang%zarot)
365 
366    ! Initialize and compute data for LDA+U
367    !paw_dmft%use_dmft=dtset%usedmft
368    !if (dtset%usepawu>0.or.dtset%useexexch>0) then
369    !  call pawpuxinit(dtset%dmatpuopt,dtset%exchmix,dtset%f4of2_sla,dtset%f6of2_sla,&
370    !    dtset%jpawu,dtset%lexexch,dtset%lpawu,cryst%ntypat,pawang,dtset%pawprtvol,&
371    !    Pawrad,pawtab,dtset%upawu,dtset%usedmft,dtset%useexexch,dtset%usepawu)
372    !end if
373    !ABI_CHECK(paw_dmft%use_dmft==0,"DMFT not available")
374    !call destroy_sc_dmft(paw_dmft)
375 
376    if (my_rank == master) call pawtab_print(pawtab, unit=std_out)
377 
378    ! Get Pawrhoij from the header of the WFK file.
379    call pawrhoij_copy(wfk0_hdr%pawrhoij,pawrhoij)
380 
381    ! Variables/arrays related to the fine FFT grid ===
382    ABI_DT_MALLOC(pawfgrtab,(cryst%natom))
383    call pawtab_get_lsize(pawtab,l_size_atm,cryst%natom,cryst%typat)
384    cplex=1
385    call pawfgrtab_init(pawfgrtab,cplex,l_size_atm,dtset%nspden,dtset%typat)
386    ABI_FREE(l_size_atm)
387 
388    usexcnhat=maxval(pawtab(:)%usexcnhat)
389    !  * 0 if Vloc in atomic data is Vbare    (Blochl s formulation)
390    !  * 1 if Vloc in atomic data is VH(tnzc) (Kresse s formulation)
391    call wrtout(std_out,sjoin("using usexcnhat= ",itoa(usexcnhat)))
392    !
393    ! Identify parts of the rectangular grid where the density has to be calculated ===
394    !optcut=0; optgr0=dtset%pawstgylm; optgr1=0; optgr2=0; optrad=1-dtset%pawstgylm
395    !if (dtset%xclevel==2 .and. usexcnhat>0) optgr1=dtset%pawstgylm
396    optcut=1; optgr0=1; optgr1=1; optgr2=1; optrad=1
397 
398    call nhatgrid(cryst%atindx1,cryst%gmet,cryst%natom,cryst%natom,cryst%nattyp,ngfftf,cryst%ntypat,&
399 &   optcut,optgr0,optgr1,optgr2,optrad,pawfgrtab,pawtab,cryst%rprimd,cryst%typat,cryst%ucvol,cryst%xred)
400 
401    call pawfgrtab_print(pawfgrtab,cryst%natom,unit=std_out,prtvol=dtset%pawprtvol)
402 
403    !ABI_MALLOC(ks_nhat,(nfftf,dtset%nspden))
404    !ks_nhat=zero
405  else
406    ABI_DT_MALLOC(pawfgrtab,(0))
407  end if !End of PAW Initialization
408 
409  select case (dtset%wfk_task)
410 
411  case (WFK_TASK_FULLBZ)
412    ! Read wfk0_path and build WFK in full BZ.
413    if (my_rank == master) then
414 
415      wfkfull_path = dtfil%fnameabo_wfk
416      if (dtset%iomode == IO_MODE_ETSF) wfkfull_path = nctk_ncify(wfkfull_path)
417      call wfk_tofullbz(wfk0_path, dtset, psps, pawtab, wfkfull_path)
418 
419      ! Write tetrahedron tables.
420      tetra = tetra_from_kptrlatt(cryst, dtset%kptopt, dtset%kptrlatt, dtset%nshiftk, &
421      dtset%shiftk, dtset%nkpt, dtset%kptns, msg, ierr)
422      if (ierr == 0) then
423        call tetra_write(tetra, dtset%nkpt, dtset%kptns, strcat(dtfil%filnam_ds(4), "_TETRA"))
424      else
425        MSG_WARNING(sjoin("Cannot produce TETRA file", ch10, msg))
426      end if
427 
428      call destroy_tetra(tetra)
429    end if
430    call xmpi_barrier(comm)
431 
432  !case ("pjdos")
433 
434  case (WFK_TASK_DDK)
435     ! calculate the DDK matrix elements using a WFK file
436     call eph_ddk(wfk0_path,dtfil,dtset,psps,pawtab,dtset%inclvkb,ngfftc,mpi_enreg,comm) 
437 
438  case (WFK_TASK_CLASSIFY)
439    ! Band classification.
440    call read_wfd()
441 
442    ABI_DT_MALLOC(esymm,(wfd%nkibz,wfd%nsppol))
443    use_paw_aeur=.False. ! should pass ngfftf but the dense mesh is not forced to be symmetric
444 
445    do spin=1,wfd%nsppol
446      do ik_ibz=1,wfd%nkibz
447        first_band = 1
448        last_band  = wfd%nband(ik_ibz,spin)
449        call classify_bands(wfd,use_paw_aeur,first_band,last_band,ik_ibz,spin,wfd%ngfft,&
450        cryst,ebands,pawtab,pawrad,pawang,psps,dtset%tolsym,esymm(ik_ibz,spin))
451      end do
452    end do
453 
454    call esymm_free(esymm)
455    ABI_DT_FREE(esymm)
456 
457  !case (WFK_TASK_UR)
458  !  ! plot KSS wavefunctions. Change bks_mask to select particular states.
459  !  ABI_MALLOC(bks_mask,(Wfd%mband,Wfd%nkibz,Wfd%nsppol))
460  !  bks_mask=.False.; bks_mask(1:4,1,1)=.True.
461  !  call wfd_plot_ur(Wfd,Cryst,Psps,Pawtab,Pawrad,ngfftf,bks_mask)
462  !  ABI_FREE(bks_mask)
463 
464  case (WFK_TASK_PAW_AEPSI)
465    ! Compute AE PAW wavefunction in real space on the dense FFT mesh.
466    call read_wfd()
467 
468    ABI_CHECK(wfd%usepaw == 1, "Not a PAW run")
469    ABI_DT_MALLOC(paw_onsite, (cryst%natom))
470    call paw_pwaves_lmn_init(paw_onsite,cryst%natom,cryst%natom,cryst%ntypat,&
471    cryst%rprimd,cryst%xcart,pawtab,pawrad,pawfgrtab)
472 
473    ! Use dense FFT mesh
474    call wfd_change_ngfft(wfd,cryst,psps,ngfftf)
475    band = 1; spin = 1; ik_ibz = 1
476 
477    ABI_MALLOC(ur_ae, (wfd%nfft*wfd%nspinor))
478    call wfd_paw_get_aeur(wfd,band,ik_ibz,spin,cryst,paw_onsite,psps,pawtab,pawfgrtab,ur_ae)
479    ABI_FREE(ur_ae)
480 
481    call paw_pwaves_lmn_free(paw_onsite)
482    ABI_DT_FREE(paw_onsite)
483 
484  !case ("paw_aeden")
485  !case ("outqmc")
486 
487  case default
488    MSG_ERROR(sjoin("Wrong task:",itoa(dtset%wfk_task)))
489  end select
490 
491  !=====================
492  !==== Free memory ====
493  !=====================
494  call crystal_free(cryst)
495  call ebands_free(ebands)
496  call wfd_free(wfd)
497  call pawfgr_destroy(pawfgr)
498  call destroy_mpi_enreg(mpi_enreg)
499  call hdr_free(wfk0_hdr)
500 
501  ! Deallocation for PAW.
502  if (dtset%usepaw==1) then
503    call pawrhoij_free(pawrhoij)
504    ABI_DT_FREE(pawrhoij)
505    call pawfgrtab_free(pawfgrtab)
506    !call paw_ij_free(paw_ij)
507    !ABI_DT_FREE(paw_ij)
508    !call paw_an_free(paw_an)
509    !ABI_DT_FREE(paw_an)
510  end if
511  ABI_DT_FREE(pawfgrtab)
512 
513  DBG_EXIT('COLL')
514 
515  contains

wfk_analyze/read_wfd [ Functions ]

[ Top ] [ wfk_analyze ] [ Functions ]

NAME

  read_wfd

FUNCTION

  Initialize the wavefunction descriptor

PARENTS

      wfk_analyze

CHILDREN

      wfd_init,wfd_read_wfk,wfd_test_ortho

SOURCE

533 subroutine read_wfd()
534 
535 
536 !This section has been created automatically by the script Abilint (TD).
537 !Do not modify the following lines by hand.
538 #undef ABI_FUNC
539 #define ABI_FUNC 'read_wfd'
540 !End of the abilint section
541 
542  implicit none
543 
544 !Arguments ------------------------------------
545 
546 !Local variables-------------------------------
547 
548 ! *************************************************************************
549 
550    ABI_MALLOC(keep_ur, (ebands%mband, ebands%nkpt, ebands%nsppol))
551    ABI_MALLOC(bks_mask, (ebands%mband, ebands%nkpt, ebands%nsppol))
552    keep_ur = .False.; bks_mask = .True.
553 
554    call wfd_init(wfd,cryst,pawtab,psps,keep_ur,dtset%paral_kgb,dummy_npw,&
555    ebands%mband,ebands%nband,ebands%nkpt,dtset%nsppol,bks_mask,&
556    dtset%nspden,dtset%nspinor,dtset%ecutsm,dtset%dilatmx,wfk0_hdr%istwfk,ebands%kptns,ngfftc,&
557    dummy_gvec,dtset%nloalg,dtset%prtvol,dtset%pawprtvol,comm,opt_ecut=ecut_eff)
558 
559    ABI_FREE(keep_ur)
560    ABI_FREE(bks_mask)
561 
562    call wfd_read_wfk(wfd,wfk0_path,IO_MODE_MPI)
563    call wfd_test_ortho(wfd,cryst,pawtab)
564 
565  end subroutine read_wfd
566 
567 end subroutine wfk_analyze