TABLE OF CONTENTS


ABINIT/m_wvl_wfsinp [ Modules ]

[ Top ] [ Modules ]

NAME

  m_wvl_wfsinp

FUNCTION

  Routines to initialize (wavelet) wavefunctions

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

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_wvl_wfsinp
28 
29  use defs_basis
30  use defs_abitypes
31  use defs_datatypes
32  use defs_wvltypes
33  use m_wffile
34  use m_abicore
35  use m_errors
36  use m_xmpi
37 
38 
39  use m_geometry,  only : xred2xcart
40  use m_abi2big,   only : wvl_occ_abi2big, wvl_occopt_abi2big, wvl_setngfft, wvl_setboxgeometry
41  use m_psolver,   only : psolver_kernel
42  use m_wvl_rwwf,  only : wvl_read
43  use m_mklocl_realspace, only : mklocl_wavelets
44  use m_wvl_wfs,          only : wvl_wfs_set, wvl_wfs_free, wvl_wfs_lr_copy
45  use m_wvl_denspot,      only : wvl_denspot_set, wvl_denspot_free
46  use m_wvl_projectors,   only : wvl_projectors_set, wvl_projectors_free
47 
48  implicit none
49 
50  private

ABINIT/wvl_wfsinp_disk [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_wfsinp_disk

FUNCTION

 This method allocates and initialises wavefunctions with values from disk.
 See wvl_wfsinp_scratch() or wvl_wfsinp_reformat() from other initialisation
 routines.

 When initialised from scratch or from disk, wvl%wfs%[h]psi comes unallocated
 and will be allocated inside this routine.
 When initialised from memory (reformating), wvl%wfs%[h]psi will be reallocated.

INPUTS

  dtset <type(dataset_type)>=input variables.
  hdr0 <type(hdr_type)>=the header of wf, den and pot files (read from restart)
  hdr <type(hdr_type)>=the header of wf, den and pot files
  mpi_enreg=informations about MPI parallelization
  option=1 for reading a file following ABINIT format, -1 for a BigDFT format.
  rprimd(3,3)=dimensional primitive translations in real space (bohr)
  wff <type(wffile_type)>= structure with informations on wf file.
  xred(3,natom)=reduced dimensionless atomic coordinates.

OUTPUT

SIDE EFFECTS

  wfs <type(wvl_projector_type)>=wavefunctions informations for wavelets.

PARENTS

      inwffil

CHILDREN

      first_orthon,wrtout,wvl_occ_abi2big,wvl_read

SOURCE

 97 subroutine wvl_wfsinp_disk(dtset, hdr0, hdr, mpi_enreg, occ, option, rprimd, wff, wfs, wvl, xred)
 98 
 99 #if defined HAVE_BIGDFT
100  use BigDFT_API, only : first_orthon,sumrho,communicate_density,plot_density
101  use dynamic_memory
102 #endif
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 'wvl_wfsinp_disk'
108 !End of the abilint section
109 
110   implicit none
111 
112 !Arguments -------------------------------
113   !scalars
114   integer, intent(in)                       :: option
115   type(dataset_type), intent(in)            :: dtset
116   type(hdr_type), intent(in)                :: hdr0
117   type(hdr_type), intent(in)                :: hdr
118   type(MPI_type), intent(in)                :: mpi_enreg
119   type(wffile_type), intent(in)             :: wff
120   type(wvl_wf_type), intent(inout)          :: wfs
121   type(wvl_internal_type), intent(inout)       :: wvl
122   !type(wvl_denspot_type), intent(inout)       :: wvl_den
123   !arrays
124   real(dp), intent(inout) :: occ(dtset%mband*dtset%nkpt*dtset%nsppol)
125   real(dp), intent(in)                      :: rprimd(3, 3)
126   real(dp), intent(in)                      :: xred(3, dtset%natom)
127 
128 !Local variables-------------------------------
129 #if defined HAVE_BIGDFT
130  logical :: wvlbigdft=.false.
131  integer :: comm,me,nproc
132  character(len = 500)  :: message
133  !real(dp), allocatable :: xcart(:,:)
134 #if defined DEBUG_MODE
135  !integer, parameter :: ndebug = 5  !5 will not work for wavelets compiling with debug=naughty
136  integer,parameter :: ndebug = 0
137 #else
138  integer, parameter :: ndebug = 0
139 #endif
140 #endif
141 
142 ! *********************************************************************
143 
144 #if defined HAVE_BIGDFT
145 
146  write(message, '(a,a)' ) ch10,' wvl_wfsinp_disk: wavefunction initialisation.'
147  call wrtout(std_out,message,'COLL')
148 
149 !If usewvl: wvlbigdft indicates that the BigDFT workflow will be followed
150  wvlbigdft=(dtset%usewvl==1.and.dtset%wvl_bigdft_comp==1)
151 
152  comm=mpi_enreg%comm_wvl
153  me=xmpi_comm_rank(comm)
154  nproc=xmpi_comm_size(comm)
155 !We allocate psi.
156 !ABI_ALLOCATE(wfs%ks%psi,( max(wfs%ks%orbs%npsidim_comp,wfs%ks%orbs%npsidim_orbs)+ndebug) )
157  wfs%ks%psi=f_malloc_ptr(max(wfs%ks%orbs%npsidim_comp,wfs%ks%orbs%npsidim_orbs)+ndebug,id='psi')
158 
159  write(message, '(a,a,a,a,I0)' ) ch10, &
160 & ' wvl_wfsinp_disk: allocate wavefunctions,', ch10, &
161 & '  size of the compressed array per proc: ', &
162 & product(shape(wfs%ks%psi))
163  call wrtout(std_out,message,'COLL')
164 
165  call wvl_read(dtset, hdr0, hdr, mpi_enreg, option, rprimd, wff, wfs, wvl, xred)
166 
167 !We orthogonalise,only for NC.
168  if(wvl%paw%usepaw==0 .and. wvlbigdft) then
169    call first_orthon(me, nproc, wfs%ks%orbs, wfs%ks%lzd, wfs%ks%comms, &
170 &   wfs%ks%psi, wfs%ks%hpsi, wfs%ks%psit, wfs%ks%orthpar,wvl%paw)
171  else
172 !  ABI_ALLOCATE(wfs%ks%hpsi,(max(wfs%ks%orbs%npsidim_orbs,wfs%ks%orbs%npsidim_comp)))
173    wfs%ks%hpsi=f_malloc_ptr(max(wfs%ks%orbs%npsidim_orbs,wfs%ks%orbs%npsidim_comp),id='hpsi')
174    if(wvl%paw%usepaw==1) then
175      ABI_ALLOCATE(wvl%paw%spsi,(max(wfs%ks%orbs%npsidim_orbs,wfs%ks%orbs%npsidim_comp)))
176    end if
177 
178 !  Set orbs%eval=-0.5.
179 !  This will be done in LDiagHam
180 !  For the moment we skip this, since hpsi is not yet calculated
181 !  and it an input argument in LDiagHam.
182    wfs%ks%orbs%eval(:)=-0.5d0
183 
184 !  Copy occupations from BigDFT objects to ABINIT
185    call wvl_occ_abi2big(dtset%mband,dtset%nkpt,dtset%nsppol,occ,2,wfs)
186 
187 
188  end if
189 
190 #else
191  BIGDFT_NOTENABLED_ERROR()
192  if (.false.) write(std_out,*) option,dtset%nstep,hdr0%ecut,hdr%ecut,mpi_enreg%me,wff%me,&
193 & wfs%ks,wvl%h(1),occ(1),rprimd(1,1),xred(1,1)
194 #endif
195 
196 !for testing
197 !!Plot the density
198 !call sumrho(wvl_den%denspot%dpbox,wfs%ks%orbs,&
199 !& wfs%GPU,& wvl%atoms%sym,&
200 !& wvl_den%denspot%rhod,wfs%psi,wvl_den%denspot%rho_psi)
201 !call communicate_density(wvl_den%denspot%dpbox,wfs%ks%orbs%nspin,&
202 !& wvl_den%denspot%rhod,wvl_den%denspot%dpcom%nscatterarr,&
203 !& wvl_den%denspot%rho_psi,wvl_den%denspot%rhov)
204 !call plot_density('electronic_density',&
205 !& me,nproc,wfs%Lzd%Glr%d%n1,wfs%Lzd%Glr%d%n2,wfs%Lzd%Glr%d%n3,&
206 !& wfs%Lzd%Glr%d%n1i,wfs%Lzd%Glr%d%n2i,wfs%Lzd%Glr%d%n3i,&
207 !& wvl_den%denspot%dpcom%nscatterarr(me,2),  &
208 !& wfs%orbs%nspin,&
209 !& wvl_den%denspot%hgrids(1),wvl_den%denspot%hgrids(2),wvl_den%denspot%hgrids(3),&
210 !& wvl%atoms,xcart,wvl_den%denspot%dpcom%ngatherarr,&
211 !& wvl_den%denspot%rhov(1+wvl_den%denspot%dpcom%nscatterarr(me,4)*wfs%Lzd%Glr%d%n1i*wfs%Lzd%Glr%d%n2i))
212 !ABI_DEALLOCATE(xcart)
213 !end of debug
214 
215 
216 end subroutine wvl_wfsinp_disk

ABINIT/wvl_wfsinp_reformat [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_wfsinp_reformat

FUNCTION

 This method allocates and initialises wavefunctions with values from disk.
 See wvl_wfsinp_scratch() or wvl_wfsinp_reformat() from other initialisation
 routines.

 When initialised from scratch or from disk, wvl%wfs%ks%[h]psi comes unallocated
 and will be allocated inside this routine.
 When initialised from memory (reformating), wvl%wfs%ks%[h]psi will be reallocated.
 The projectors are also recomputed.

 The scalar arrays should be reallocated using dtset%nfft after a call to
 this routine.

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

      mover

CHILDREN

      copy_old_wavefunctions,deallocate_wfd,first_orthon
      local_potential_dimensions,nullify_gaussian_basis,psolver_kernel
      reformatmywaves,wrtout,wvl_denspot_free,wvl_denspot_set
      wvl_projectors_free,wvl_projectors_set,wvl_setboxgeometry,wvl_setngfft
      wvl_wfs_free,wvl_wfs_lr_copy,wvl_wfs_set,xred2xcart

SOURCE

257 subroutine wvl_wfsinp_reformat(dtset, mpi_enreg, psps, rprimd, wvl, xred, xred_old)
258 
259 #if defined HAVE_BIGDFT
260  use BigDFT_API, only : copy_old_wavefunctions, reformatmywaves, first_orthon, &
261 & deallocate_wfd, wavefunctions_descriptors, deallocate_lr, &
262 & local_potential_dimensions, copy_coulomb_operator, &
263 & deallocate_coulomb_operator, nullify_gaussian_basis
264  use dynamic_memory
265 #endif
266 
267 !This section has been created automatically by the script Abilint (TD).
268 !Do not modify the following lines by hand.
269 #undef ABI_FUNC
270 #define ABI_FUNC 'wvl_wfsinp_reformat'
271 !End of the abilint section
272 
273   implicit none
274 
275 !Arguments ------------------------------------
276   type(dataset_type), intent(inout)      :: dtset
277   type(MPI_type), intent(inout)          :: mpi_enreg
278   type(pseudopotential_type), intent(in) :: psps
279   type(wvl_data), intent(inout)          :: wvl
280   real(dp), intent(inout)                :: rprimd(3,3)
281   real(dp), intent(inout)                :: xred_old(3, dtset%natom)
282   real(dp), intent(inout)                :: xred(3, dtset%natom)
283 
284 !Local variables-------------------------------
285 #if defined HAVE_BIGDFT
286   integer                  :: itypat
287   integer                  :: nSize_old(3)
288   real(dp)                 :: hgrid_old(3)
289   real(dp), allocatable    :: xcart(:,:), xcart_old(:,:)
290   real(dp), pointer        :: psi_old(:), eigen_old(:)
291   integer :: comm,me,nproc,icoulomb
292   type(coulomb_operator)::kernel
293   type(wavefunctions_descriptors) :: keys_old
294   character(len=500)       :: message
295 #if defined DEBUG_MODE
296  !integer, parameter :: ndebug = 5  !5 will not work for wavelets compiling with debug=naughty
297  integer,parameter :: ndebug = 0
298 #else
299  integer, parameter :: ndebug = 0
300 #endif
301 #endif
302 
303 ! *********************************************************************
304 
305 #if defined HAVE_BIGDFT
306 
307  write(message, '(a,a)' ) ch10,&
308 & ' wvl_wfsinp_reformat: reformat the wavefunctions.'
309  call wrtout(std_out, message, 'COLL')
310 
311  comm=mpi_enreg%comm_wvl
312  me=xmpi_comm_rank(comm)
313  nproc=xmpi_comm_size(comm)
314 !Convert input xred_old (reduced coordinates) to xcart_old (cartesian)
315  ABI_ALLOCATE(xcart_old,(3, dtset%natom))
316  call xred2xcart(dtset%natom, rprimd, xcart_old, xred_old)
317 
318 !Copy current to old.
319  ABI_ALLOCATE(eigen_old,(wvl%wfs%ks%orbs%norb))
320  eigen_old = wvl%wfs%ks%orbs%eval
321  hgrid_old = wvl%descr%h
322  call copy_old_wavefunctions(nproc, wvl%wfs%ks%orbs, &
323 & wvl%descr%Glr%d%n1, wvl%descr%Glr%d%n2, wvl%descr%Glr%d%n3, &
324 & wvl%wfs%ks%lzd%Glr%wfd, wvl%wfs%ks%psi, nSize_old(1), nSize_old(2), nSize_old(3), &
325 & keys_old, psi_old)
326 !Patch because copy_old_wavefunctions() free wvl%wfs%ks%lzd%Glr%wfd but don't nullify it.
327  nullify(wvl%wfs%ks%lzd%glr%wfd%keyglob)
328  nullify(wvl%wfs%ks%lzd%glr%wfd%keygloc)
329  nullify(wvl%wfs%ks%lzd%glr%wfd%keyvloc)
330  nullify(wvl%wfs%ks%lzd%glr%wfd%keyvglob)
331 
332 !We deallocate the previous projectors.
333  call wvl_projectors_free(wvl%projectors)
334 
335 !Deallocate old wavefunctions
336  call wvl_wfs_free(wvl%wfs)
337 
338 !Deallocate old denspot
339  call wvl_denspot_free(wvl%den)
340 
341 !We change the box geometry.
342  call wvl_setBoxGeometry(dtset%prtvol, psps%gth_params%radii_cf, rprimd, xred, &
343 & wvl%descr, dtset%wvl_crmult, dtset%wvl_frmult)
344  call wvl_denspot_set(wvl%den, psps%gth_params, dtset%ixc, dtset%natom, dtset%nsppol, &
345 & rprimd, wvl%descr, dtset%wvl_crmult, dtset%wvl_frmult, comm, xred)
346  if (wvl%descr%atoms%astruct%geocode == "F") then
347    icoulomb = 1
348  else if (wvl%descr%atoms%astruct%geocode == "S") then
349    icoulomb = 2
350  else
351    icoulomb = 0
352  end if
353 !calculation of the Poisson kernel anticipated to reduce memory peak for small systems
354  call psolver_kernel( wvl%den%denspot%dpbox%hgrids, 2, icoulomb, me, kernel, &
355 & comm, wvl%den%denspot%dpbox%ndims, nproc, dtset%nscforder)
356  ! Shallow copy of the kernel (still owned by ABINIT).
357  wvl%den%denspot%pkernel = kernel
358  wvl%den%denspot%pkernelseq = kernel
359 !Associate the denspot distribution into mpi_enreg.
360  mpi_enreg%nscatterarr => wvl%den%denspot%dpbox%nscatterarr
361  mpi_enreg%ngatherarr => wvl%den%denspot%dpbox%ngatherarr
362  mpi_enreg%ngfft3_ionic = wvl%den%denspot%dpbox%n3pi
363  call wvl_setngfft(me, dtset%mgfft, dtset%nfft, &
364 & dtset%ngfft, nproc, wvl%den%denspot%dpbox%ndims(1), &
365 & wvl%den%denspot%dpbox%ndims(2), &
366 & wvl%den%denspot%dpbox%ndims(3),wvl%den%denspot%dpbox%n3d)
367 
368 !We copy the geometry structure.
369  call wvl_wfs_lr_copy(wvl%wfs, wvl%descr)
370 !Reallocate them with new size.
371  call wvl_wfs_set(dtset%strprecon,dtset%spinmagntarget, dtset%kpt, me, dtset%natom, sum(dtset%nband), dtset%nkpt, &
372 & nproc, dtset%nspinor, dtset%nsppol, dtset%nwfshist, dtset%occ_orig(:,1), psps, rprimd, &
373 & wvl%wfs, dtset%wtk, wvl%descr, dtset%wvl_crmult, dtset%wvl_frmult, xred)
374 
375 !Recopy old eval for precond.
376  wvl%wfs%ks%orbs%eval = eigen_old
377  ABI_DEALLOCATE(eigen_old)
378 
379 !We allocate psi.
380 !ABI_ALLOCATE(wvl%wfs%ks%psi,( max(wvl%wfs%ks%orbs%npsidim_comp,wvl%wfs%ks%orbs%npsidim_orbs)+ndebug) )
381  wvl%wfs%ks%psi=f_malloc_ptr(max(wvl%wfs%ks%orbs%npsidim_comp,wvl%wfs%ks%orbs%npsidim_orbs)+ndebug,id='psi')
382  write(message, '(a,a,a,a,I0)' ) ch10, &
383 & ' wvl_wfsinp_reformat: allocate wavefunctions,', ch10, &
384 & '  size of the compressed array per proc: ', &
385 & product(shape(wvl%wfs%ks%psi))
386  call wrtout(std_out,message,'COLL')
387 
388 !Convert input xred (reduced coordinates) to xcart (cartesian)
389  ABI_ALLOCATE(xcart,(3, dtset%natom))
390  call xred2xcart(dtset%natom, rprimd, xcart, xred)
391 
392 !We transfer the old wavefunctions to the new ones.
393  call reformatmywaves(me, wvl%wfs%ks%orbs, wvl%descr%atoms, &
394 & hgrid_old(1), hgrid_old(2), hgrid_old(3), nSize_old(1), nSize_old(2), &
395 & nSize_old(3), xcart_old, keys_old, psi_old, wvl%descr%h(1), wvl%descr%h(2), &
396 & wvl%descr%h(3), wvl%descr%Glr%d%n1, wvl%descr%Glr%d%n2, wvl%descr%Glr%d%n3, xcart, &
397 & wvl%wfs%ks%lzd%Glr%wfd, wvl%wfs%ks%psi)
398  ABI_DEALLOCATE(xcart)
399  ABI_DEALLOCATE(xcart_old)
400 
401 !We free the old descriptors and arrays.
402  ABI_DEALLOCATE(psi_old)
403  call deallocate_wfd(keys_old)
404 
405  call local_potential_dimensions(me,wvl%wfs%ks%lzd,wvl%wfs%ks%orbs,wvl%den%denspot%xc,&
406 & wvl%den%denspot%dpbox%ngatherarr(0,1))
407 
408 !it seems that the table "wvl%projectors%G" is no more used
409 !but it's not allocated -> fortran runtime error
410 #if defined HAVE_BIGDFT
411  ABI_DATATYPE_ALLOCATE(wvl%projectors%G,(dtset%ntypat))
412  do itypat=1,dtset%ntypat
413    call nullify_gaussian_basis(wvl%projectors%G(itypat))
414  end do
415 #endif
416 
417 !Reallocate projectors for the new positions.
418  call wvl_projectors_set(me, dtset%natom, wvl%projectors, psps, rprimd, &
419 & wvl%wfs, wvl%descr, dtset%wvl_frmult, xred)
420 
421 !Orthogonilise new wavefunctions.
422  call first_orthon(me, nproc, wvl%wfs%ks%orbs, wvl%wfs%ks%lzd, wvl%wfs%ks%comms, &
423 & wvl%wfs%ks%psi, wvl%wfs%ks%hpsi, wvl%wfs%ks%psit, wvl%wfs%ks%orthpar,wvl%descr%paw)
424 
425 #else
426  BIGDFT_NOTENABLED_ERROR()
427  if (.false.) write(std_out,*) dtset%nstep,mpi_enreg%me,psps%npsp,wvl%wfs%ks,rprimd(1,1),&
428 & xred_old(1,1),xred(1,1)
429 #endif
430 
431 end subroutine wvl_wfsinp_reformat

ABINIT/wvl_wfsinp_scratch [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_wfsinp_scratch

FUNCTION

 This method allocates and initialises wavefunctions with values from input guess.
 See wvl_wfsinp_disk() or wvl_wfsinp_reformat() from other initialisation
 routines.

 When initialised from scratch or from disk, wvl%wfs%[h]psi comes unallocated
 and will be allocated inside this routine.
 When initialised from memory (reformating), wvl%wfs%[h]psi will be reallocated.

INPUTS

  dtset <type(dataset_type)>=input variables.
  hdr0 <type(hdr_type)>=the header of wf, den and pot files (read from restart)
  hdr <type(hdr_type)>=the header of wf, den and pot files
  ireadwf=1 for reading from file, 0 otherwise.
  mpi_enreg=informations about MPI parallelization
  option=1 for reading a file following ABINIT format, -1 for a BigDFT format.
  rprimd(3,3)=dimensional primitive translations in real space (bohr)
  wff <type(wffile_type)>= structure with informations on wf file.
  xred(3,natom)=reduced dimensionless atomic coordinates

OUTPUT

SIDE EFFECTS

  wvl <type(wvl_data)>=wavefunctions & projectors informations for wavelets.

PARENTS

      inwffil

CHILDREN

      input_wf_diag,mklocl_wavelets,wrtout,wvl_occ_abi2big,wvl_occopt_abi2big
      xred2xcart

SOURCE

472 subroutine wvl_wfsinp_scratch(dtset, mpi_enreg, occ, rprimd, wvl, xred)
473 
474 #if defined HAVE_BIGDFT
475  use BigDFT_API, only : createIonicPotential, input_wf_diag, gaussian_basis, &
476       & input_variables, calculate_rhocore, deallocate_Lzd_except_Glr, INPUT_IG_OFF,&
477       & SMEARING_DIST_ERF, PSPCODE_PAW
478 #endif
479 
480 !This section has been created automatically by the script Abilint (TD).
481 !Do not modify the following lines by hand.
482 #undef ABI_FUNC
483 #define ABI_FUNC 'wvl_wfsinp_scratch'
484 !End of the abilint section
485 
486   implicit none
487 
488 !Arguments -------------------------------
489   !scalars
490   type(dataset_type), intent(in)        :: dtset
491   type(MPI_type), intent(inout)         :: mpi_enreg
492   type(wvl_data), intent(inout)         :: wvl
493   !arrays
494   real(dp), intent(inout) :: occ(dtset%mband*dtset%nkpt*dtset%nsppol)
495   real(dp), intent(in)                  :: rprimd(3, 3)
496   real(dp), intent(in)                  :: xred(3, dtset%natom)
497 
498 !Local variables-------------------------------
499 #if defined HAVE_BIGDFT
500   character(len = 500)  :: message
501   integer               :: comm,me,nproc
502   integer               :: iscf_local
503   integer               :: nvirt
504   integer               :: ii,shift_vpsp,size_vpsp
505   logical               :: onlywf=.false. ! find the wavefunctions and return
506   logical               :: wvlbigdft=.false.
507   real(dp), allocatable :: xcart(:,:)
508   real(dp), allocatable :: rhor(:,:)
509   real(dp), pointer     :: vpsp(:)
510   real(dp):: elecfield(3)
511   type(gaussian_basis) :: Gvirt
512   type(input_variables) :: in  ! To be removed, waiting for BigDFT upgrade
513 #endif
514 
515 ! *********************************************************************
516 
517 #if defined HAVE_BIGDFT
518 
519  elecfield=zero
520 
521 !If usewvl: wvlbigdft indicates that the BigDFT workflow will be followed
522  wvlbigdft=(dtset%usewvl==1.and.dtset%wvl_bigdft_comp==1)
523 
524  write(message, '(a,a)' ) ch10,&
525 & ' wvl_wfsinp_scratch: wavefunction initialisation.'
526  call wrtout(std_out,message,'COLL')
527 
528  comm=mpi_enreg%comm_wvl
529  me=xmpi_comm_rank(comm)
530  nproc=xmpi_comm_size(comm)
531 !Store xcart for each atom
532  ABI_ALLOCATE(xcart,(3, dtset%natom))
533  call xred2xcart(dtset%natom, rprimd, xcart, xred)
534 
535 !We allocate temporary arrays for rho and vpsp.
536 !allocate ionic potential
537  if (wvl%den%denspot%dpbox%n3pi > 0) then
538    size_vpsp=wvl%descr%Glr%d%n1i*wvl%descr%Glr%d%n2i*wvl%den%denspot%dpbox%n3pi
539    shift_vpsp=wvl%den%denspot%dpbox%ndims(1)*wvl%den%denspot%dpbox%ndims(2) &
540 &   *wvl%den%denspot%dpbox%nscatterarr(me,4)
541    ABI_ALLOCATE(vpsp,(size_vpsp+shift_vpsp))
542  else
543    ABI_ALLOCATE(vpsp,(1))
544  end if
545 
546  if(.not. onlywf) then
547    ABI_ALLOCATE(rhor,(dtset%nfft,dtset%nspden))
548    call mklocl_wavelets(elecfield, xcart, mpi_enreg, dtset%natom, &
549 &   dtset%nfft, dtset%nspden, 1, rprimd, vpsp, &
550 &   wvl%den, wvl%descr, xcart)
551    ABI_DEALLOCATE(rhor)
552  end if
553 
554 ! IMPORTANT: onlywf=.true. does not work yet, do not change this:
555 ! if(.not. wvlbigdft) onlywf=.true. !do not apply Hamiltonian inside input_wf_diag.
556  if(dtset%usepaw==1) wvl%descr%atoms%npspcode(:)=wvl%descr%npspcode_paw_init_guess(:)
557  iscf_local=dtset%iscf
558  if(.not. wvlbigdft) iscf_local=0  !important to have good occ values
559 
560 !This routine allocates psi, hpsi and psit.
561  nvirt = 0
562  in%linear = INPUT_IG_OFF
563  in%nspin = dtset%nsppol
564  in%exctxpar = wvl%descr%exctxpar
565  in%itrpmax = dtset%nnsclo
566  in%iscf = iscf_local !dtset%iscf
567  in%Tel = dtset%tsmear
568 ! if (dtset%iscf == 0) in%Tel = zero
569  if (iscf_local == 0) in%Tel = zero
570  in%SIC = wvl%wfs%ks%SIC
571  in%orthpar = wvl%wfs%ks%orthpar
572 
573  in%occopt=dtset%occopt
574  if(dtset%occopt>2 .and. dtset%occopt<7)then
575    call wvl_occopt_abi2big(in%occopt,in%occopt,1)
576  else
577 !  This will be used only for the initial wavefunctions:
578    in%occopt=SMEARING_DIST_ERF
579  end if
580 
581 !Note: check if all required in "in" is passed.
582 !remove argument "in" from input_wf_diag
583  call input_wf_diag(me, nproc, &
584 & wvl%descr%atoms, wvl%den%denspot,&
585 & wvl%wfs%ks%orbs, nvirt, wvl%wfs%ks%comms, &
586 & wvl%wfs%ks%lzd, wvl%e%energs, xcart, &
587 & wvl%projectors%nlpsp, &
588 & dtset%ixc, wvl%wfs%ks%psi, wvl%wfs%ks%hpsi, wvl%wfs%ks%psit, &
589 & Gvirt, in%nspin, wvl%wfs%GPU, in, onlywf, wvl%projectors%G, wvl%descr%paw)
590 
591 !This provisory: wvl%descr%paw could be passed as optional
592 !to input_wf_diag to allocate spsi inside this routine
593  if(dtset%usepaw==1) then
594    wvl%descr%atoms%npspcode(:)=PSPCODE_PAW
595    ABI_ALLOCATE(wvl%descr%paw%spsi,(max(wvl%wfs%ks%orbs%npsidim_orbs,wvl%wfs%ks%orbs%npsidim_comp)))
596    do ii=1,size(wvl%wfs%ks%psi)
597      wvl%descr%paw%spsi(ii)=wvl%wfs%ks%psi(ii)
598      wvl%wfs%ks%hpsi(ii)=wvl%wfs%ks%psi(ii)
599    end do
600  end if
601 
602  if(wvlbigdft ) then
603 !  Copy occupations from BigDFT objects to ABINIT
604    call wvl_occ_abi2big(dtset%mband,dtset%nkpt,dtset%nsppol,occ,2,wvl%wfs)
605  else
606 !  Copy occupations from ABINIT to BigDFT objects
607    call wvl_occ_abi2big(dtset%mband,dtset%nkpt,dtset%nsppol,occ,1,wvl%wfs)
608  end if
609 
610  write(message, '(a)' ) &
611 & '  | wavefunctions have been calculated.'
612  call wrtout(std_out,message,'COLL')
613 
614  ABI_DEALLOCATE(xcart)
615  ABI_DEALLOCATE(vpsp)
616 
617 #else
618  BIGDFT_NOTENABLED_ERROR()
619  if (.false.) write(std_out,*) dtset%nstep,mpi_enreg%me,wvl%wfs%ks,occ(1),rprimd(1,1),xred(1,1)
620 #endif
621 
622 end subroutine wvl_wfsinp_scratch