TABLE OF CONTENTS


ABINIT/m_pawfgr [ Modules ]

[ Top ] [ Modules ]

NAME

  m_pawfgr

FUNCTION

  This module contains the definition of the pawfgr_type structured datatype,
  as well as related functions and methods.
  pawfgr_type variables define Fine rectangular GRid parameters and related data.

COPYRIGHT

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

INPUTS

OUTPUT

NOTES

  * Routines tagged with "@type_name" are strongly connected to the definition of the data type. 
    Strongly connected means that the proper functioning of the implementation relies on the 
    assumption that the tagged procedure is consistent with the type declaration.
    Every time a developer changes the structure "type_name" adding new entries, he/she has to make sure 
    that all the strongly connected routines are changed accordingly to accommodate the modification of the data type
    Typical examples of strongly connected routines are creation, destruction or reset methods.

PARENTS

CHILDREN

SOURCE

35 #if defined HAVE_CONFIG_H
36 #include "config.h"
37 #endif
38 
39 #include "abi_common.h"
40 
41 MODULE m_pawfgr
42 
43  use defs_basis
44  use m_errors
45  use m_abicore
46  use m_xmpi
47 
48  use defs_abitypes, only : dataset_type
49  use m_kg,       only : getcut
50 
51  implicit none
52 
53  private
54 
55 !public procedures.
56  public :: pawfgr_init
57  public :: pawfgr_destroy
58  public :: pawfgr_nullify
59  public :: indgrid

m_paw_ij/pawfgr_nullify [ Functions ]

[ Top ] [ m_paw_ij ] [ Functions ]

NAME

  pawfgr_nullify

FUNCTION

  Nullify pointers and flags in a pawfgr structure

SIDE EFFECTS

  Pawfgr<type(pawfgr_type)>=Fine GRid parameters and related data. Nullified in output

PARENTS

      m_rec

CHILDREN

SOURCE

366 subroutine pawfgr_nullify(Pawfgr)
367 
368 
369 !This section has been created automatically by the script Abilint (TD).
370 !Do not modify the following lines by hand.
371 #undef ABI_FUNC
372 #define ABI_FUNC 'pawfgr_nullify'
373 !End of the abilint section
374 
375  implicit none
376 
377 !Arguments ------------------------------------
378 !arrays
379  type(Pawfgr_type),intent(inout) :: Pawfgr
380 
381 !Local variables-------------------------------
382 
383 ! *************************************************************************
384 
385  !@Pawfgr_type
386 
387   nullify(Pawfgr%coatofin)
388   nullify(Pawfgr%fintocoa)
389 
390   Pawfgr%usefinegrid=0
391 
392 end subroutine pawfgr_nullify

m_pawfgr/indgrid [ Functions ]

[ Top ] [ m_pawfgr ] [ Functions ]

NAME

 indgrid

FUNCTION

 Calculate the correspondance between the coarse grid and the fine grid

COPYRIGHT

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

INPUTS

 nfftc=total number of FFt grid=n1*n2*n3 for the coarse grid
 nfftf=total number of FFt grid=n1*n2*n3 for the fine grid
 ngfftc(18)=contain all needed information about 3D FFT, for the coarse grid,
        see ~abinit/doc/variables/vargs.htm#ngfft
 ngfftf(18)=contain all needed information about 3D FFT, for the fine grid,
        see ~abinit/doc/variables/vargs.htm#ngfft

OUTPUT

 coatofin(nfftc)= index of the points of the coarse grid on the fine grid
 fintocoa(nfftf)=index of the points of the fine grid on the
   coarse grid (=0 if the point of the fine grid does not belong to
   the coarse grid).

PARENTS

      fourier_interpol,m_pawfgr,m_rec

CHILDREN

SOURCE

433 subroutine indgrid(coatofin,fintocoa,nfftc,nfftf,ngfftc,ngfftf)
434 
435 
436 !This section has been created automatically by the script Abilint (TD).
437 !Do not modify the following lines by hand.
438 #undef ABI_FUNC
439 #define ABI_FUNC 'indgrid'
440 !End of the abilint section
441 
442  implicit none
443 
444 !Arguments ------------------------------------
445 !scalars
446  integer,intent(in) :: nfftc,nfftf
447 !arrays
448  integer,intent(in) :: ngfftc(18),ngfftf(18)
449  integer,intent(out) :: coatofin(nfftc),fintocoa(nfftf)
450 
451 !Local variables-------------------------------
452 !scalars
453  integer :: i1,i2,i3,if1,if2,if3,ii,ing,n1c,n1f,n2c,n2f,n3c,n3f,narg1,narg2
454 !arrays
455  integer :: id(3)
456  integer,allocatable :: gc(:,:),gf(:,:)
457  character(len=500) :: msg
458 
459 ! *************************************************************************
460 !
461 
462  DBG_ENTER("COLL")
463 
464  n1c=ngfftc(1);n2c=ngfftc(2);n3c=ngfftc(3)
465  n1f=ngfftf(1);n2f=ngfftf(2);n3f=ngfftf(3)
466 
467  ABI_ALLOCATE(gc,(3,max(n1c,n2c,n3c)))
468  do ii=1,3
469    id(ii)=ngfftc(ii)/2+2
470    do ing=1,ngfftc(ii)
471      gc(ii,ing)=ing-(ing/id(ii))*ngfftc(ii)-1
472    end do
473  end do
474 
475  ABI_ALLOCATE(gf,(3,max(n1f,n2f,n3f)))
476  do ii=1,3
477    id(ii)=ngfftf(ii)/2+2
478    do ing=1,ngfftf(ii)
479      gf(ii,ing)=ing-(ing/id(ii))*ngfftf(ii)-1
480    end do
481  end do
482 
483  coatofin=0;fintocoa=0
484  do i1=1,n1c
485    do if1=1,n1f
486      if(gc(1,i1)==gf(1,if1)) then
487        do i2=1,n2c
488          do if2=1,n2f
489            if(gc(2,i2)==gf(2,if2)) then
490              do i3=1,n3c
491                do if3=1,n3f
492                  if(gc(3,i3)==gf(3,if3)) then
493                    narg1=i1+n1c*(i2-1+n2c*(i3-1))
494                    narg2=if1+n1f*(if2-1+n2f*(if3-1))
495                    coatofin(narg1)=narg2
496                    fintocoa(narg2)=narg1
497                    exit
498                  end if
499                end do
500              end do
501              exit ! To avoid N_fine * N_coarse scaling
502            end if
503          end do
504        end do
505        exit ! To avoid N_fine * N_coarse scaling
506      end if
507    end do
508  end do
509 
510 !Check coatofin to make sure there are no zeros!
511  do ii=1,ubound(coatofin,1)
512    if (coatofin(ii)==0) then
513      msg = 'A zero was found in coatofin. Check that the fine FFT mesh is finer in each dimension than the coarse FFT mesh.'
514      MSG_ERROR(msg)
515    end if
516  end do
517 
518  ABI_DEALLOCATE(gf)
519  ABI_DEALLOCATE(gc)
520 
521  DBG_EXIT("COLL")
522 
523 end subroutine indgrid

m_pawfgr/pawfgr_destroy [ Functions ]

[ Top ] [ m_pawfgr ] [ Functions ]

NAME

  pawfgr_destroy

FUNCTION

  Deallocate pointers and nullify flags in a pawfgr structure

SIDE EFFECTS

  pawfgr<type(pawfgr_type)>= Fine GRid parameters and related data

PARENTS

      bethe_salpeter,eph,fourier_interpol,gstate,m_rec,respfn,screening,sigma
      wfk_analyze

CHILDREN

SOURCE

310 subroutine pawfgr_destroy(Pawfgr)
311 
312 
313 !This section has been created automatically by the script Abilint (TD).
314 !Do not modify the following lines by hand.
315 #undef ABI_FUNC
316 #define ABI_FUNC 'pawfgr_destroy'
317 !End of the abilint section
318 
319  implicit none
320 
321 !Arguments ------------------------------------
322 !arrays
323  type(Pawfgr_type),intent(inout) :: Pawfgr
324 
325 !Local variables-------------------------------
326 
327 ! *************************************************************************
328 
329  DBG_ENTER("COLL")
330 
331 !@Pawfgr_type
332  
333   if (associated(Pawfgr%coatofin))  then
334     ABI_DEALLOCATE(Pawfgr%coatofin)
335   end if
336   if (associated(Pawfgr%fintocoa))  then
337     ABI_DEALLOCATE(Pawfgr%fintocoa)
338   end if
339 
340   Pawfgr%usefinegrid=0
341 
342  DBG_EXIT("COLL")
343 
344 end subroutine pawfgr_destroy

m_pawfgr/pawfgr_init [ Functions ]

[ Top ] [ m_pawfgr ] [ Functions ]

NAME

 pawfgr_init

FUNCTION

  Initialize a pawfgr_type datatype, reporting also info on the mesh
  according to the method used (norm-conserving PSP or PAW)

INPUTS

  k0(3)=input k vector for k+G sphere
  Dtset <type(dataset_type)>=all input variables for this dataset
   %dilatmx
   %usepaw
   %usewvl
   %natom
   %ngfft
   %ngfftdg
   %nfft
   %mgfft
   %mgfftdg
   %dilatmx
   %pawecutdg
   %ecut
  gmet(3,3)=reciprocal space metric (bohr^-2)

OUTPUT

  ecut_eff=effective energy cutoff (hartree) for coarse planewave basis sphere
  ecutdg_eff=effective energy cutoff (hartree) for dense planewave basis sphere
  gsqcutc_eff=(PAW) Fourier cutoff on G^2 for "large sphere" of radius double for the coarse FFT grid
  nfftf=(effective) number of FFT grid points (for this proc), for dense FFT mesh
  mgfftf=maximum size of 1D FFTs, for dense FFT mesh
  ngfftc(18),ngfftf(18)=contain all needed information about 3D FFT, for coarse and dense FFT mesh, resp.
                        see ~abinit/doc/variables/vargs.htm#ngfft
  Pawfgr<pawfgr_type>=For PAW, Fine rectangular GRid parameters and related data

PARENTS

      bethe_salpeter,eph,gstate,respfn,screening,sigma,wfk_analyze

CHILDREN

SOURCE

168 subroutine pawfgr_init(Pawfgr,Dtset,mgfftf,nfftf,ecut_eff,ecutdg_eff,ngfftc,ngfftf,&
169 &                      gsqcutc_eff,gsqcutf_eff,gmet,k0) ! optional
170 
171 
172 !This section has been created automatically by the script Abilint (TD).
173 !Do not modify the following lines by hand.
174 #undef ABI_FUNC
175 #define ABI_FUNC 'pawfgr_init'
176 !End of the abilint section
177 
178  implicit none
179 
180 !Arguments ------------------------------------
181 !scalars
182  integer,intent(out) :: nfftf,mgfftf
183  real(dp),intent(out) :: ecut_eff,ecutdg_eff
184  real(dp),intent(out),optional :: gsqcutf_eff,gsqcutc_eff
185  type(dataset_type),intent(in) :: Dtset
186  type(Pawfgr_type),intent(out) :: Pawfgr
187 !arrays
188  real(dp),intent(in),optional :: gmet(3,3)
189  integer,intent(out) :: ngfftc(18),ngfftf(18)
190  real(dp),intent(in),optional :: k0(3)
191 
192 !Local variables-------------------------------
193  integer :: ii,nfftc_tot,nfftf_tot
194  real(dp) :: boxcut,boxcutc
195  character(len=500) :: msg
196 
197 !************************************************************************
198 
199  DBG_ENTER("COLL")
200 
201  !@Pawfgr_type
202 
203  if ((present(gsqcutc_eff).or.present(gsqcutf_eff)).and.&
204 &    ((.not.present(gmet)).or.(.not.present(k0)))) then
205    msg='To compute gsqcut[c,f]_eff, both k0 and gmet must be present as argument !'
206    MSG_BUG(msg)
207  end if
208 
209  ngfftc(:)=Dtset%ngfft(:)
210 
211  SELECT CASE (Dtset%usepaw)
212 
213  CASE (0)
214   ! === Norm-conserving pseudopotentials ===
215   nfftf=Dtset%nfft ; mgfftf=Dtset%mgfft ; ngfftf(:)=Dtset%ngfft(:)
216   Pawfgr%usefinegrid=0
217   ABI_ALLOCATE(Pawfgr%coatofin,(0))
218   ABI_ALLOCATE(Pawfgr%fintocoa,(0))
219   ecut_eff  =Dtset%ecut*Dtset%dilatmx**2
220   ecutdg_eff=ecut_eff
221 
222  CASE (1)
223   ! == PAW calculation ===
224   if (any(Dtset%ngfftdg(1:3)/=Dtset%ngfft(1:3)) .and. Dtset%usewvl==0) then
225 ! if (Dtset%pawecutdg>=1.0000001_dp*Dtset%ecut .and. Dtset%usewvl==0) then
226    ! * Use fine FFT grid generated according to pawecutdg.
227    nfftf=Dtset%nfftdg ; mgfftf=Dtset%mgfftdg ; ngfftf(:)=Dtset%ngfftdg(:)
228    nfftc_tot =ngfftc(1)*ngfftc(2)*ngfftc(3)
229    nfftf_tot =ngfftf(1)*ngfftf(2)*ngfftf(3)
230    Pawfgr%usefinegrid=1
231    ABI_ALLOCATE(Pawfgr%coatofin,(nfftc_tot))
232    ABI_ALLOCATE(Pawfgr%fintocoa,(nfftf_tot))
233    call indgrid(Pawfgr%coatofin,Pawfgr%fintocoa,nfftc_tot,nfftf_tot,ngfftc,ngfftf)
234   else
235    ! * Do not use fine FFT mesh. Simple transfer that can be done in parallel with only local info.
236    nfftf=Dtset%nfft ; mgfftf=Dtset%mgfft ; ngfftf(:)=Dtset%ngfft(:)
237    Pawfgr%usefinegrid=0
238    ABI_ALLOCATE(Pawfgr%coatofin,(Dtset%nfft))
239    ABI_ALLOCATE(Pawfgr%fintocoa,(Dtset%nfft))
240    do ii=1,Dtset%nfft
241     Pawfgr%coatofin(ii)=ii ; Pawfgr%fintocoa(ii)=ii
242    end do
243   end if
244   ecutdg_eff=Dtset%pawecutdg*Dtset%dilatmx**2
245   ecut_eff  =Dtset%ecut*Dtset%dilatmx**2
246 
247  CASE DEFAULT
248   write(msg,'(a,i4)')' Wrong value of usepaw: ',Dtset%usepaw
249   MSG_BUG(msg)
250  END SELECT
251 
252 ! == Store useful dimensions in Pawfgr ===
253  Pawfgr%nfftc=Dtset%nfft ; Pawfgr%mgfftc=Dtset%mgfft ; Pawfgr%ngfftc(:)=Dtset%ngfft(:)
254  Pawfgr%nfft=nfftf       ; Pawfgr%mgfft=mgfftf       ; Pawfgr%ngfft (:)=ngfftf(:)
255 
256  !
257  ! === Get boxcut for given gmet, ngfft, and ecut (center at k0) ===
258  !     boxcut=ratio of basis sphere diameter to fft box side
259  boxcut=-one
260  if (Dtset%usepaw==1) then
261    if (present(gsqcutc_eff)) then
262      write(msg,'(2a)')ch10,' Coarse grid specifications (used for wave-functions):'
263      call wrtout(std_out,msg,'COLL')
264      call getcut(boxcutc,ecut_eff,gmet,gsqcutc_eff,Dtset%iboxcut,std_out,k0,ngfftc)
265    end if
266    if (present(gsqcutf_eff)) then
267      write(msg,'(2a)')ch10,' Fine grid specifications (used for densities):'
268      call wrtout(std_out,msg,'COLL')
269      call getcut(boxcut,ecutdg_eff,gmet,gsqcutf_eff,Dtset%iboxcut,std_out,k0,ngfftf)
270    end if
271  else if (present(gsqcutc_eff)) then
272    call getcut(boxcut,ecut_eff,gmet,gsqcutc_eff,Dtset%iboxcut,std_out,k0,ngfftc)
273    gsqcutf_eff=gsqcutc_eff
274  end if
275  !
276  ! === Check that boxcut>=2 if intxc=1; otherwise intxc must be set=0 ===
277  if (boxcut>=zero .and. boxcut<two .and. Dtset%intxc==1) then
278    write(msg,'(a,es12.4,5a)')&
279 &   ' boxcut=',boxcut,' is < 2.0  => intxc must be 0;',ch10,&
280 &   ' Need larger ngfft to use intxc=1.',ch10,&
281 &   ' Action : you could increase ngfft, or decrease ecut, or put intxc=0.'
282    MSG_ERROR(msg)
283  end if
284 
285  DBG_EXIT("COLL")
286 
287 end subroutine pawfgr_init

m_pawfgr/pawfgr_type [ Types ]

[ Top ] [ m_pawfgr ] [ Types ]

NAME

 pawfgr_type

FUNCTION

 For PAW, Fine GRid parameters and related data

SOURCE

 74  type,public :: pawfgr_type
 75 
 76 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines,
 77 ! declared in another part of ABINIT, that might need to take into account your modification.
 78 
 79 !Integer scalars
 80 
 81   integer :: mgfft, nfft
 82    ! Values of mffft and nfft for the fine grid:
 83    !   mgfft= max(ngfft(i)) [max. size of 1D FFT grid]
 84    !   nfft=ngfft1*ngfft2*ngfft3 [number of pts in the FFT box]
 85 
 86   integer :: mgfftc, nfftc
 87    ! Values of mffft and nfft for the COARSE grid:
 88    !   mgfftc= max(ngfftc(i)) [max. size of 1D FFT grid]
 89    !   nfftc=ngfftc1*ngfftc2*ngfftc3 [number of pts in the FFT box]
 90 
 91   integer :: usefinegrid
 92    ! Flag: =1 if a double-grid is used to convert spherical data
 93    !       to Fourier grid. =0 otherwise
 94 
 95 !Integer arrays
 96 
 97   ! MGTODO: Replace with allocatable
 98   integer, pointer :: coatofin(:)
 99    ! coatofin(nfftc)
100    ! Index of the points of the coarse grid on the fine grid
101 
102   integer, pointer :: fintocoa(:)
103    ! fintocoa(nfft)
104    ! Index of the points of the fine grid on the coarse grid
105    !  (=0 if the point of the fine grid does not belong to the coarse grid)
106 
107   integer :: ngfft(18)
108    ! ngfft(1:18)=integer array with FFT box dimensions and other
109    ! information on FFTs, for the fine rectangular grid
110 
111   integer :: ngfftc(18)
112    ! ngfft(1:18)=integer array with FFT box dimensions and other
113    ! information on FFTs, for the COARSE rectangular grid
114 
115  end type pawfgr_type