TABLE OF CONTENTS


ABINIT/derf_ab [ Functions ]

[ Top ] [ Functions ]

NAME

 derf_ab

FUNCTION

 Some wrappers for BigDFT which uses different names for the same routines.

INPUTS

OUTPUT

SIDE EFFECTS

PARENTS

      mklocl_realspace,mklocl_wavelets

CHILDREN

SOURCE

315 subroutine derf_ab(derf_yy,yy)
316 
317  use m_special_funcs,  only : abi_derf
318 
319 !This section has been created automatically by the script Abilint (TD).
320 !Do not modify the following lines by hand.
321 #undef ABI_FUNC
322 #define ABI_FUNC 'derf_ab'
323 !End of the abilint section
324 
325  implicit none
326 
327 !Arguments ------------------------------------
328 !scalars
329  real(dp),intent(in) :: yy
330  real(dp),intent(out) :: derf_yy
331 
332 !Local variables-------------------------------
333 
334 ! *********************************************************************
335 
336  derf_yy = abi_derf(yy)
337 
338 end subroutine derf_ab

ABINIT/derfcf [ Functions ]

[ Top ] [ Functions ]

NAME

 derfcf

FUNCTION

 Some wrappers for BigDFT which uses different names for the same routines.

INPUTS

OUTPUT

SIDE EFFECTS

PARENTS

CHILDREN

SOURCE

271 subroutine derfcf(derfc_yy,yy)
272 
273  use m_special_funcs,  only : abi_derfc
274 
275 !This section has been created automatically by the script Abilint (TD).
276 !Do not modify the following lines by hand.
277 #undef ABI_FUNC
278 #define ABI_FUNC 'derfcf'
279 !End of the abilint section
280 
281  implicit none
282 !Arguments ------------------------------------
283 !scalars
284  real(dp),intent(in) :: yy
285  real(dp),intent(out) :: derfc_yy
286 !Local variables-------------------------------
287 
288 ! *********************************************************************
289 
290  derfc_yy = abi_derfc(yy)
291 
292 end subroutine derfcf

ABINIT/m_wvl_wfs [ Modules ]

[ Top ] [ Modules ]

NAME

 m_wvl_wfs

FUNCTION

COPYRIGHT

  Copyright (C) 1998-2018 ABINIT group (DC)
  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

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

ABINIT/wvl_wfs_free [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_wfs_free

FUNCTION

 Freeing routine.

INPUTS

OUTPUT

SIDE EFFECTS

  wfs <type(wvl_wf_type)>=wavefunctions informations in a wavelet basis.

PARENTS

      gstate,wvl_wfsinp_reformat

CHILDREN

      deallocate_comms,deallocate_lr,deallocate_lzd_except_glr
      deallocate_orbs,f_free_ptr

SOURCE

367 subroutine wvl_wfs_free(wfs)
368 
369  use defs_wvltypes
370 #if defined HAVE_BIGDFT
371  use BigDFT_API, only: deallocate_Lzd_except_Glr, deallocate_lr, &
372       & deallocate_orbs, deallocate_comms
373  use dynamic_memory
374 #endif
375 
376 !This section has been created automatically by the script Abilint (TD).
377 !Do not modify the following lines by hand.
378 #undef ABI_FUNC
379 #define ABI_FUNC 'wvl_wfs_free'
380 !End of the abilint section
381 
382  implicit none
383 
384 !Arguments ------------------------------------
385 !scalars
386  type(wvl_wf_type),intent(inout) :: wfs
387 
388 !Local variables -------------------------
389 
390 ! *********************************************************************
391 
392 #if defined HAVE_BIGDFT
393  call deallocate_Lzd_except_Glr(wfs%ks%lzd)
394  call deallocate_lr(wfs%ks%lzd%Glr)
395  call deallocate_orbs(wfs%ks%orbs)
396  call deallocate_comms(wfs%ks%comms)
397  if (associated(wfs%ks%orbs%eval))  then
398    ABI_DEALLOCATE(wfs%ks%orbs%eval)
399  end if
400  ABI_DATATYPE_DEALLOCATE(wfs%ks%confdatarr)
401 
402  if (associated(wfs%ks%psi)) then
403    call f_free_ptr(wfs%ks%psi)
404  end if
405  if (associated(wfs%ks%hpsi)) then
406    call f_free_ptr(wfs%ks%hpsi)
407  end if
408  if (associated(wfs%ks%psit)) then
409    call f_free_ptr(wfs%ks%psit)
410  end if
411 
412 #else
413  BIGDFT_NOTENABLED_ERROR()
414  if (.false.) write(std_out,*) wfs%ks
415 #endif
416 
417 end subroutine wvl_wfs_free

ABINIT/wvl_wfs_set [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_wfs_set

FUNCTION

 Compute the access keys for the wavefunctions when the positions
 of the atoms are given.

 For memory occupation optimisation reasons, the wavefunctions are not allocated
 here. See the initialisation routines wvl_wfsinp_disk(), wvl_wfsinp_scratch()
 and wvl_wfsinp_reformat() to do it. After allocation, use wvl_wfs_free()
 to deallocate all stuff (descriptors and arrays).

INPUTS

  dtset <type(dataset_type)>=internal variables used by wavelets, describing
   | wvl_internal=desciption of the wavelet box.
   | natom=number of atoms.
  mpi_enreg=informations about MPI parallelization
  psps <type(pseudopotential_type)>=variables related to pseudopotentials
  rprimd(3,3)=dimensional primitive translations in real space (bohr)
  xred(3,natom)=reduced dimensionless atomic coordinates

OUTPUT

  wfs <type(wvl_projector_type)>=wavefunctions informations for wavelets.
   | keys=its access keys for compact storage.

PARENTS

      gstate,wvl_wfsinp_reformat

CHILDREN

SOURCE

 79 subroutine wvl_wfs_set(alphadiis, spinmagntarget, kpt, me, natom, nband, nkpt, nproc, nspinor, &
 80 &  nsppol, nwfshist, occ, psps, rprimd, wfs, wtk, wvl, wvl_crmult, wvl_frmult, xred)
 81 
 82  use defs_datatypes
 83  use defs_wvltypes
 84 
 85  use m_geometry, only : xred2xcart
 86 #if defined HAVE_BIGDFT
 87  use BigDFT_API, only: createWavefunctionsDescriptors, orbitals_descriptors, &
 88        & orbitals_communicators, allocate_diis_objects, &
 89        & input_variables, check_linear_and_create_Lzd, check_communications, &
 90        & INPUT_IG_OFF, nullify_locreg_descriptors
 91 #endif
 92 
 93 !This section has been created automatically by the script Abilint (TD).
 94 !Do not modify the following lines by hand.
 95 #undef ABI_FUNC
 96 #define ABI_FUNC 'wvl_wfs_set'
 97 !End of the abilint section
 98 
 99  implicit none
100 
101 !Arguments ------------------------------------
102 !scalars
103  integer, intent(in) :: natom, nkpt, nsppol, nspinor, nband, nwfshist,me,nproc
104  real(dp), intent(in) :: spinmagntarget, wvl_crmult, wvl_frmult, alphadiis
105  type(pseudopotential_type),intent(in) :: psps
106  type(wvl_wf_type),intent(out) :: wfs
107  type(wvl_internal_type), intent(in) :: wvl
108 !arrays
109  real(dp), intent(in) :: kpt(3,nkpt)
110  real(dp), intent(in) :: wtk(nkpt), occ(:)
111  real(dp),intent(in) :: rprimd(3,3),xred(3,natom)
112 
113 !Local variables-------------------------------
114 #if defined HAVE_BIGDFT
115 !scalars
116  integer :: ii,idata, norb, norbu, norbd
117  logical :: parallel
118  character(len=500) :: message
119  type(input_variables) :: in  ! To be removed, waiting for BigDFT upgrade
120 !arrays
121  real(dp), allocatable :: kpt_(:,:)
122  real(dp),allocatable :: xcart(:,:)
123 #endif
124 
125 ! *********************************************************************
126 
127 #if defined HAVE_BIGDFT
128 
129  parallel = (nproc > 1)
130 
131 !Consistency checks, are all pseudo true GTH pseudo with geometric informations?
132 !Skip for PAW case: we do not have GTH parameters
133  do idata = 1, psps%npsp, 1
134    if (.not. psps%gth_params%set(idata) .and. psps%usepaw==0) then
135      write(message, '(a,a,a,a,I0,a,a,a)' ) ch10,&
136 &     ' wvl_wfs_set:  consistency checks failed,', ch10, &
137 &     '  no GTH parameters found for type number ', idata, '.', ch10, &
138 &     '  Check your input pseudo files.'
139      MSG_ERROR(message)
140    end if
141    if (.not. psps%gth_params%hasGeometry(idata)) then
142      write(message, '(a,a,a,a,a,a)' ) ch10,&
143 &     ' wvl_wfs_set:  consistency checks failed,', ch10, &
144 &     '  the given GTH parameters has no geometry informations.', ch10, &
145 &     '  Upgrade your input pseudo files to GTH with geometric informations.'
146      MSG_ERROR(message)
147    end if
148  end do
149 
150 !Store xcart for each atom
151  ABI_ALLOCATE(xcart,(3, natom))
152  call xred2xcart(natom, rprimd, xcart, xred)
153 
154 !Nullify possibly unset pointers
155  nullify(wfs%ks%psi)
156  nullify(wfs%ks%hpsi)
157  nullify(wfs%ks%psit)
158 
159 !Static allocations.
160  norb = nband / nkpt
161  norbu = 0
162  norbd = 0
163  if (nsppol == 2) then
164    if (spinmagntarget < -real(90, dp)) then
165      norbu = min(norb / 2, norb)
166    else
167      norbu = min(norb / 2 + int(spinmagntarget), norb)
168    end if
169    norbd = norb - norbu
170  else
171    norbu = norb
172    norbd = 0
173  end if
174  ABI_ALLOCATE(kpt_, (3, nkpt))
175  do ii = 1, nkpt
176    kpt_(:,ii) = kpt(:,ii) / (/ rprimd(1,1), rprimd(2,2), rprimd(3,3) /) * two_pi
177  end do
178 
179  call orbitals_descriptors(me, nproc,norb,norbu,norbd,nsppol,nspinor, &
180 & nkpt,kpt_,wtk,wfs%ks%orbs,.false.)
181  ABI_DEALLOCATE(kpt_)
182 !We copy occ_orig to wfs%ks%orbs%occup
183  wfs%ks%orbs%occup(1:norb * nkpt) = occ(1:norb * nkpt)
184 !We allocate the eigen values storage.
185  ABI_ALLOCATE(wfs%ks%orbs%eval,(wfs%ks%orbs%norb * wfs%ks%orbs%nkpts))
186 
187  write(message, '(a,a)' ) ch10,&
188 & ' wvl_wfs_set: Create access keys for wavefunctions.'
189  call wrtout(std_out,message,'COLL')
190 
191  call nullify_locreg_descriptors(wfs%ks%lzd%Glr)
192  wfs%ks%lzd%Glr = wvl%Glr
193  call createWavefunctionsDescriptors(me, wvl%h(1), wvl%h(2), wvl%h(3), &
194 & wvl%atoms, xcart, psps%gth_params%radii_cf, &
195 & wvl_crmult, wvl_frmult, wfs%ks%lzd%Glr)
196 !The memory is not allocated there for memory occupation optimisation reasons.
197 
198  call orbitals_communicators(me,nproc,wfs%ks%lzd%Glr,wfs%ks%orbs,wfs%ks%comms)
199 
200  write(message, '(a,2I8)' ) &
201 & '  | all orbitals have coarse segments, elements:', &
202 & wfs%ks%lzd%Glr%wfd%nseg_c, wfs%ks%lzd%Glr%wfd%nvctr_c
203  call wrtout(std_out,message,'COLL')
204  write(message, '(a,2I8)' ) &
205 & '  | all orbitals have fine   segments, elements:', &
206 & wfs%ks%lzd%Glr%wfd%nseg_f, 7 * wfs%ks%lzd%Glr%wfd%nvctr_f
207  call wrtout(std_out,message,'COLL')
208 
209 !allocate arrays necessary for DIIS convergence acceleration
210  call allocate_diis_objects(nwfshist,alphadiis,&
211 & sum(wfs%ks%comms%ncntt(0:nproc-1)), wfs%ks%orbs%nkptsp, wfs%ks%orbs%nspinor, &
212 & wfs%ks%diis)
213 
214  ABI_DATATYPE_ALLOCATE(wfs%ks%confdatarr, (wfs%ks%orbs%norbp))
215  call default_confinement_data(wfs%ks%confdatarr,wfs%ks%orbs%norbp)
216 
217  call check_linear_and_create_Lzd(me,nproc,INPUT_IG_OFF,wfs%ks%lzd,&
218 & wvl%atoms,wfs%ks%orbs,nsppol,xcart)
219  wfs%ks%lzd%hgrids = wvl%h
220 
221 !check the communication distribution
222  call check_communications(me,nproc,wfs%ks%orbs,wfs%ks%Lzd,wfs%ks%comms)
223 
224 !Deallocations
225  ABI_DEALLOCATE(xcart)
226 
227 !DEBUG
228  write(std_out,*) 'wvl_wfs_set: TODO, update BigDFT sic_input_variables_default()'
229 !ENDDEBUG
230  call sic_input_variables_default(in)
231 
232  wfs%ks%SIC                  = in%SIC
233  wfs%ks%exctxpar             = "OP2P"
234  wfs%ks%c_obj                = 0
235  wfs%ks%orthpar%directDiag   = .true.
236  wfs%ks%orthpar%norbpInguess = 5
237  wfs%ks%orthpar%bsLow        = 300
238  wfs%ks%orthpar%bsUp         = 800
239  wfs%ks%orthpar%methOrtho    = 0
240  wfs%ks%orthpar%iguessTol    = 1.d-4
241 
242 #else
243  BIGDFT_NOTENABLED_ERROR()
244  if (.false.) write(std_out,*) natom,nkpt,nsppol,nspinor,nband,nwfshist,me,nproc,&
245 & spinmagntarget,wvl_crmult,wvl_frmult,alphadiis,psps%npsp,wfs%ks,wvl%h(1),&
246 & kpt(1,1),wtk(1),occ(1),rprimd(1,1),xred(1,1)
247 #endif
248 
249 end subroutine wvl_wfs_set

defs_wvltypes/wvl_wfs_lr_copy [ Functions ]

[ Top ] [ defs_wvltypes ] [ Functions ]

NAME

 wvl_wfs_lr_copy

FUNCTION

 Copy the wvl%Glr datastructure geometry part to wfs%Glr.

INPUTS

 wvl <type(wvl_internal_type)> = input localisation region

OUTPUT

 wfs <type(wvl_wf_type)> = output localistaion region

PARENTS

      gstate,wvl_wfsinp_reformat

CHILDREN

SOURCE

442 subroutine wvl_wfs_lr_copy(wfs, wvl)
443 
444  use defs_wvltypes
445 
446 !This section has been created automatically by the script Abilint (TD).
447 !Do not modify the following lines by hand.
448 #undef ABI_FUNC
449 #define ABI_FUNC 'wvl_wfs_lr_copy'
450 !End of the abilint section
451 
452  implicit none
453 
454 !Arguments ------------------------------------
455 !scalars
456   type(wvl_internal_type), intent(in)  :: wvl
457   type(wvl_wf_type), intent(inout)     :: wfs
458 !arrays
459 
460 !Local variables-------------------------------
461 
462 ! *********************************************************************
463 
464 #if defined HAVE_BIGDFT
465 !Use global localization region for the moment.
466  wfs%ks%lzd%Glr%geocode    = wvl%Glr%geocode
467  wfs%ks%lzd%Glr%hybrid_on  = wvl%Glr%hybrid_on
468  wfs%ks%lzd%Glr%ns1        = wvl%Glr%ns1
469  wfs%ks%lzd%Glr%ns2        = wvl%Glr%ns2
470  wfs%ks%lzd%Glr%ns3        = wvl%Glr%ns3
471  wfs%ks%lzd%Glr%nsi1       = wvl%Glr%nsi1
472  wfs%ks%lzd%Glr%nsi2       = wvl%Glr%nsi2
473  wfs%ks%lzd%Glr%nsi3       = wvl%Glr%nsi3
474  wfs%ks%lzd%Glr%d%n1       = wvl%Glr%d%n1
475  wfs%ks%lzd%Glr%d%n2       = wvl%Glr%d%n2
476  wfs%ks%lzd%Glr%d%n3       = wvl%Glr%d%n3
477  wfs%ks%lzd%Glr%d%nfl1     = wvl%Glr%d%nfl1
478  wfs%ks%lzd%Glr%d%nfu1     = wvl%Glr%d%nfu1
479  wfs%ks%lzd%Glr%d%nfl2     = wvl%Glr%d%nfl2
480  wfs%ks%lzd%Glr%d%nfu2     = wvl%Glr%d%nfu2
481  wfs%ks%lzd%Glr%d%nfl3     = wvl%Glr%d%nfl3
482  wfs%ks%lzd%Glr%d%nfu3     = wvl%Glr%d%nfu3
483  wfs%ks%lzd%Glr%d%n1i      = wvl%Glr%d%n1i
484  wfs%ks%lzd%Glr%d%n2i      = wvl%Glr%d%n2i
485  wfs%ks%lzd%Glr%d%n3i      = wvl%Glr%d%n3i
486  wfs%ks%lzd%Glr%outofzone  = wvl%Glr%outofzone
487 
488 #else
489  BIGDFT_NOTENABLED_ERROR()
490  if (.false.) write(std_out,*) wvl%h(1),wfs%ks
491 #endif
492 
493 end subroutine wvl_wfs_lr_copy