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

SOURCE

283 subroutine derf_ab(derf_yy,yy)
284 
285  use m_special_funcs,  only : abi_derf
286  implicit none
287 
288 !Arguments ------------------------------------
289 !scalars
290  real(dp),intent(in) :: yy
291  real(dp),intent(out) :: derf_yy
292 
293 !Local variables-------------------------------
294 
295 ! *********************************************************************
296 
297  derf_yy = abi_derf(yy)
298 
299 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

SOURCE

251 subroutine derfcf(derfc_yy,yy)
252 
253  use m_special_funcs,  only : abi_derfc
254  implicit none
255 !Arguments ------------------------------------
256 !scalars
257  real(dp),intent(in) :: yy
258  real(dp),intent(out) :: derfc_yy
259 !Local variables-------------------------------
260 
261 ! *********************************************************************
262 
263  derfc_yy = abi_derfc(yy)
264 
265 end subroutine derfcf

ABINIT/m_wvl_wfs [ Modules ]

[ Top ] [ Modules ]

NAME

 m_wvl_wfs

FUNCTION

COPYRIGHT

  Copyright (C) 1998-2022 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 .

SOURCE

15 #if defined HAVE_CONFIG_H
16 #include "config.h"
17 #endif
18 
19 #include "abi_common.h"
20 
21 module m_wvl_wfs
22 
23  use defs_basis
24  use m_errors
25  use m_abicore
26 
27  use defs_datatypes, only : pseudopotential_type
28 
29  implicit none
30 
31  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 information in a wavelet basis.

SOURCE

321 subroutine wvl_wfs_free(wfs)
322 
323  use defs_wvltypes
324 #if defined HAVE_BIGDFT
325  use BigDFT_API, only: deallocate_Lzd_except_Glr, deallocate_lr, &
326       & deallocate_orbs, deallocate_comms
327  use dynamic_memory
328 #endif
329  implicit none
330 
331 !Arguments ------------------------------------
332 !scalars
333  type(wvl_wf_type),intent(inout) :: wfs
334 
335 !Local variables -------------------------
336 
337 ! *********************************************************************
338 
339 #if defined HAVE_BIGDFT
340  call deallocate_Lzd_except_Glr(wfs%ks%lzd)
341  call deallocate_lr(wfs%ks%lzd%Glr)
342  call deallocate_orbs(wfs%ks%orbs)
343  call deallocate_comms(wfs%ks%comms)
344  if (associated(wfs%ks%orbs%eval))  then
345    ABI_FREE(wfs%ks%orbs%eval)
346  end if
347  ABI_FREE(wfs%ks%confdatarr)
348 
349  if (associated(wfs%ks%psi)) then
350    call f_free_ptr(wfs%ks%psi)
351  end if
352  if (associated(wfs%ks%hpsi)) then
353    call f_free_ptr(wfs%ks%hpsi)
354  end if
355  if (associated(wfs%ks%psit)) then
356    call f_free_ptr(wfs%ks%psit)
357  end if
358 
359 #else
360  BIGDFT_NOTENABLED_ERROR()
361  if (.false.) write(std_out,*) wfs%ks
362 #endif
363 
364 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=information 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 information for wavelets.
   | keys=its access keys for compact storage.

SOURCE

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

SOURCE

384 subroutine wvl_wfs_lr_copy(wfs, wvl)
385 
386  use defs_wvltypes
387  implicit none
388 
389 !Arguments ------------------------------------
390 !scalars
391   type(wvl_internal_type), intent(in)  :: wvl
392   type(wvl_wf_type), intent(inout)     :: wfs
393 !arrays
394 
395 !Local variables-------------------------------
396 
397 ! *********************************************************************
398 
399 #if defined HAVE_BIGDFT
400 !Use global localization region for the moment.
401  wfs%ks%lzd%Glr%geocode    = wvl%Glr%geocode
402  wfs%ks%lzd%Glr%hybrid_on  = wvl%Glr%hybrid_on
403  wfs%ks%lzd%Glr%ns1        = wvl%Glr%ns1
404  wfs%ks%lzd%Glr%ns2        = wvl%Glr%ns2
405  wfs%ks%lzd%Glr%ns3        = wvl%Glr%ns3
406  wfs%ks%lzd%Glr%nsi1       = wvl%Glr%nsi1
407  wfs%ks%lzd%Glr%nsi2       = wvl%Glr%nsi2
408  wfs%ks%lzd%Glr%nsi3       = wvl%Glr%nsi3
409  wfs%ks%lzd%Glr%d%n1       = wvl%Glr%d%n1
410  wfs%ks%lzd%Glr%d%n2       = wvl%Glr%d%n2
411  wfs%ks%lzd%Glr%d%n3       = wvl%Glr%d%n3
412  wfs%ks%lzd%Glr%d%nfl1     = wvl%Glr%d%nfl1
413  wfs%ks%lzd%Glr%d%nfu1     = wvl%Glr%d%nfu1
414  wfs%ks%lzd%Glr%d%nfl2     = wvl%Glr%d%nfl2
415  wfs%ks%lzd%Glr%d%nfu2     = wvl%Glr%d%nfu2
416  wfs%ks%lzd%Glr%d%nfl3     = wvl%Glr%d%nfl3
417  wfs%ks%lzd%Glr%d%nfu3     = wvl%Glr%d%nfu3
418  wfs%ks%lzd%Glr%d%n1i      = wvl%Glr%d%n1i
419  wfs%ks%lzd%Glr%d%n2i      = wvl%Glr%d%n2i
420  wfs%ks%lzd%Glr%d%n3i      = wvl%Glr%d%n3i
421  wfs%ks%lzd%Glr%outofzone  = wvl%Glr%outofzone
422 
423 #else
424  BIGDFT_NOTENABLED_ERROR()
425  if (.false.) write(std_out,*) wvl%h(1),wfs%ks
426 #endif
427 
428 end subroutine wvl_wfs_lr_copy