TABLE OF CONTENTS
ABINIT/m_wvl_wfsinp [ Modules ]
NAME
m_wvl_wfsinp
FUNCTION
Routines to initialize (wavelet) wavefunctions
COPYRIGHT
Copyright (C) 1998-2024 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
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_wvl_wfsinp 23 24 use defs_basis 25 use defs_wvltypes 26 use m_wffile 27 use m_abicore 28 use m_errors 29 use m_xmpi 30 use m_hdr 31 use m_dtset 32 33 use defs_datatypes, only : pseudopotential_type 34 use defs_abitypes, only : MPI_type 35 use m_geometry, only : xred2xcart 36 use m_abi2big, only : wvl_occ_abi2big, wvl_occopt_abi2big, wvl_setngfft, wvl_setboxgeometry 37 use m_psolver, only : psolver_kernel 38 use m_wvl_rwwf, only : wvl_read 39 use m_mklocl_realspace, only : mklocl_wavelets 40 use m_wvl_wfs, only : wvl_wfs_set, wvl_wfs_free, wvl_wfs_lr_copy 41 use m_wvl_denspot, only : wvl_denspot_set, wvl_denspot_free 42 use m_wvl_projectors, only : wvl_projectors_set, wvl_projectors_free 43 44 implicit none 45 46 private
ABINIT/wvl_wfsinp_disk [ 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=information 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 information on wf file. xred(3,natom)=reduced dimensionless atomic coordinates.
OUTPUT
SIDE EFFECTS
wfs <type(wvl_projector_type)>=wavefunctions information for wavelets.
SOURCE
87 subroutine wvl_wfsinp_disk(dtset, hdr0, hdr, mpi_enreg, occ, option, rprimd, wff, wfs, wvl, xred) 88 89 #if defined HAVE_BIGDFT 90 use BigDFT_API, only : first_orthon,sumrho,communicate_density,plot_density 91 use dynamic_memory 92 #endif 93 94 !Arguments ------------------------------- 95 !scalars 96 integer, intent(in) :: option 97 type(dataset_type), intent(in) :: dtset 98 type(hdr_type), intent(in) :: hdr0 99 type(hdr_type), intent(in) :: hdr 100 type(MPI_type), intent(in) :: mpi_enreg 101 type(wffile_type), intent(in) :: wff 102 type(wvl_wf_type), intent(inout) :: wfs 103 type(wvl_internal_type), intent(inout) :: wvl 104 !type(wvl_denspot_type), intent(inout) :: wvl_den 105 !arrays 106 real(dp), intent(inout) :: occ(dtset%mband*dtset%nkpt*dtset%nsppol) 107 real(dp), intent(in) :: rprimd(3, 3) 108 real(dp), intent(in) :: xred(3, dtset%natom) 109 110 !Local variables------------------------------- 111 #if defined HAVE_BIGDFT 112 logical :: wvlbigdft=.false. 113 integer :: comm,me,nproc 114 character(len = 500) :: message 115 !real(dp), allocatable :: xcart(:,:) 116 #if defined DEBUG_MODE 117 !integer, parameter :: ndebug = 5 !5 will not work for wavelets compiling with debug=naughty 118 integer,parameter :: ndebug = 0 119 #else 120 integer, parameter :: ndebug = 0 121 #endif 122 #endif 123 124 ! ********************************************************************* 125 126 #if defined HAVE_BIGDFT 127 128 write(message, '(a,a)' ) ch10,' wvl_wfsinp_disk: wavefunction initialisation.' 129 call wrtout(std_out,message,'COLL') 130 131 !If usewvl: wvlbigdft indicates that the BigDFT workflow will be followed 132 wvlbigdft=(dtset%usewvl==1.and.dtset%wvl_bigdft_comp==1) 133 134 comm=mpi_enreg%comm_wvl 135 me=xmpi_comm_rank(comm) 136 nproc=xmpi_comm_size(comm) 137 !We allocate psi. 138 !ABI_MALLOC(wfs%ks%psi,( max(wfs%ks%orbs%npsidim_comp,wfs%ks%orbs%npsidim_orbs)+ndebug) ) 139 wfs%ks%psi=f_malloc_ptr(max(wfs%ks%orbs%npsidim_comp,wfs%ks%orbs%npsidim_orbs)+ndebug,id='psi') 140 141 write(message, '(a,a,a,a,I0)' ) ch10, & 142 & ' wvl_wfsinp_disk: allocate wavefunctions,', ch10, & 143 & ' size of the compressed array per proc: ', & 144 & product(shape(wfs%ks%psi)) 145 call wrtout(std_out,message,'COLL') 146 147 call wvl_read(dtset, hdr0, hdr, mpi_enreg, option, rprimd, wff, wfs, wvl, xred) 148 149 !We orthogonalise,only for NC. 150 if(wvl%paw%usepaw==0 .and. wvlbigdft) then 151 call first_orthon(me, nproc, wfs%ks%orbs, wfs%ks%lzd, wfs%ks%comms, & 152 & wfs%ks%psi, wfs%ks%hpsi, wfs%ks%psit, wfs%ks%orthpar,wvl%paw) 153 else 154 ! ABI_MALLOC(wfs%ks%hpsi,(max(wfs%ks%orbs%npsidim_orbs,wfs%ks%orbs%npsidim_comp))) 155 wfs%ks%hpsi=f_malloc_ptr(max(wfs%ks%orbs%npsidim_orbs,wfs%ks%orbs%npsidim_comp),id='hpsi') 156 if(wvl%paw%usepaw==1) then 157 ABI_MALLOC(wvl%paw%spsi,(max(wfs%ks%orbs%npsidim_orbs,wfs%ks%orbs%npsidim_comp))) 158 end if 159 160 ! Set orbs%eval=-0.5. 161 ! This will be done in LDiagHam 162 ! For the moment we skip this, since hpsi is not yet calculated 163 ! and it an input argument in LDiagHam. 164 wfs%ks%orbs%eval(:)=-0.5d0 165 166 ! Copy occupations from BigDFT objects to ABINIT 167 call wvl_occ_abi2big(dtset%mband,dtset%nkpt,dtset%nsppol,occ,2,wfs) 168 169 170 end if 171 172 #else 173 BIGDFT_NOTENABLED_ERROR() 174 if (.false.) write(std_out,*) option,dtset%nstep,hdr0%ecut,hdr%ecut,mpi_enreg%me,wff%me,& 175 & wfs%ks,wvl%h(1),occ(1),rprimd(1,1),xred(1,1) 176 #endif 177 178 !for testing 179 !!Plot the density 180 !call sumrho(wvl_den%denspot%dpbox,wfs%ks%orbs,& 181 !& wfs%GPU,& wvl%atoms%sym,& 182 !& wvl_den%denspot%rhod,wfs%psi,wvl_den%denspot%rho_psi) 183 !call communicate_density(wvl_den%denspot%dpbox,wfs%ks%orbs%nspin,& 184 !& wvl_den%denspot%rhod,wvl_den%denspot%dpcom%nscatterarr,& 185 !& wvl_den%denspot%rho_psi,wvl_den%denspot%rhov) 186 !call plot_density('electronic_density',& 187 !& me,nproc,wfs%Lzd%Glr%d%n1,wfs%Lzd%Glr%d%n2,wfs%Lzd%Glr%d%n3,& 188 !& wfs%Lzd%Glr%d%n1i,wfs%Lzd%Glr%d%n2i,wfs%Lzd%Glr%d%n3i,& 189 !& wvl_den%denspot%dpcom%nscatterarr(me,2), & 190 !& wfs%orbs%nspin,& 191 !& wvl_den%denspot%hgrids(1),wvl_den%denspot%hgrids(2),wvl_den%denspot%hgrids(3),& 192 !& wvl%atoms,xcart,wvl_den%denspot%dpcom%ngatherarr,& 193 !& wvl_den%denspot%rhov(1+wvl_den%denspot%dpcom%nscatterarr(me,4)*wfs%Lzd%Glr%d%n1i*wfs%Lzd%Glr%d%n2i)) 194 !ABI_FREE(xcart) 195 !end of debug 196 197 198 end subroutine wvl_wfsinp_disk
ABINIT/wvl_wfsinp_reformat [ 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
SOURCE
229 subroutine wvl_wfsinp_reformat(dtset, mpi_enreg, psps, rprimd, wvl, xred, xred_old) 230 231 #if defined HAVE_BIGDFT 232 use BigDFT_API, only : copy_old_wavefunctions, reformatmywaves, first_orthon, & 233 & deallocate_wfd, wavefunctions_descriptors, deallocate_lr, & 234 & local_potential_dimensions, copy_coulomb_operator, & 235 & deallocate_coulomb_operator, nullify_gaussian_basis 236 use dynamic_memory 237 #endif 238 239 !Arguments ------------------------------------ 240 type(dataset_type), intent(inout) :: dtset 241 type(MPI_type), intent(inout) :: mpi_enreg 242 type(pseudopotential_type), intent(in) :: psps 243 type(wvl_data), intent(inout) :: wvl 244 real(dp), intent(inout) :: rprimd(3,3) 245 real(dp), intent(inout) :: xred_old(3, dtset%natom) 246 real(dp), intent(inout) :: xred(3, dtset%natom) 247 248 !Local variables------------------------------- 249 #if defined HAVE_BIGDFT 250 integer :: itypat 251 integer :: nSize_old(3) 252 real(dp) :: hgrid_old(3) 253 real(dp), allocatable :: xcart(:,:), xcart_old(:,:) 254 real(dp), pointer :: psi_old(:), eigen_old(:) 255 integer :: comm,me,nproc,icoulomb 256 type(coulomb_operator)::kernel 257 type(wavefunctions_descriptors) :: keys_old 258 character(len=500) :: message 259 #if defined DEBUG_MODE 260 !integer, parameter :: ndebug = 5 !5 will not work for wavelets compiling with debug=naughty 261 integer,parameter :: ndebug = 0 262 #else 263 integer, parameter :: ndebug = 0 264 #endif 265 #endif 266 267 ! ********************************************************************* 268 269 #if defined HAVE_BIGDFT 270 271 write(message, '(a,a)' ) ch10,& 272 & ' wvl_wfsinp_reformat: reformat the wavefunctions.' 273 call wrtout(std_out, message, 'COLL') 274 275 comm=mpi_enreg%comm_wvl 276 me=xmpi_comm_rank(comm) 277 nproc=xmpi_comm_size(comm) 278 !Convert input xred_old (reduced coordinates) to xcart_old (cartesian) 279 ABI_MALLOC(xcart_old,(3, dtset%natom)) 280 call xred2xcart(dtset%natom, rprimd, xcart_old, xred_old) 281 282 !Copy current to old. 283 ABI_MALLOC(eigen_old,(wvl%wfs%ks%orbs%norb)) 284 eigen_old = wvl%wfs%ks%orbs%eval 285 hgrid_old = wvl%descr%h 286 call copy_old_wavefunctions(nproc, wvl%wfs%ks%orbs, & 287 & wvl%descr%Glr%d%n1, wvl%descr%Glr%d%n2, wvl%descr%Glr%d%n3, & 288 & wvl%wfs%ks%lzd%Glr%wfd, wvl%wfs%ks%psi, nSize_old(1), nSize_old(2), nSize_old(3), & 289 & keys_old, psi_old) 290 !Patch because copy_old_wavefunctions() free wvl%wfs%ks%lzd%Glr%wfd but don't nullify it. 291 nullify(wvl%wfs%ks%lzd%glr%wfd%keyglob) 292 nullify(wvl%wfs%ks%lzd%glr%wfd%keygloc) 293 nullify(wvl%wfs%ks%lzd%glr%wfd%keyvloc) 294 nullify(wvl%wfs%ks%lzd%glr%wfd%keyvglob) 295 296 !We deallocate the previous projectors. 297 call wvl_projectors_free(wvl%projectors) 298 299 !Deallocate old wavefunctions 300 call wvl_wfs_free(wvl%wfs) 301 302 !Deallocate old denspot 303 call wvl_denspot_free(wvl%den) 304 305 !We change the box geometry. 306 call wvl_setBoxGeometry(dtset%prtvol, psps%gth_params%radii_cf, rprimd, xred, & 307 & wvl%descr, dtset%wvl_crmult, dtset%wvl_frmult) 308 call wvl_denspot_set(wvl%den, psps%gth_params, dtset%ixc, dtset%natom, dtset%nsppol, & 309 & rprimd, wvl%descr, dtset%wvl_crmult, dtset%wvl_frmult, comm, xred) 310 if (wvl%descr%atoms%astruct%geocode == "F") then 311 icoulomb = 1 312 else if (wvl%descr%atoms%astruct%geocode == "S") then 313 icoulomb = 2 314 else 315 icoulomb = 0 316 end if 317 !calculation of the Poisson kernel anticipated to reduce memory peak for small systems 318 call psolver_kernel( wvl%den%denspot%dpbox%hgrids, 2, icoulomb, me, kernel, & 319 & comm, wvl%den%denspot%dpbox%ndims, nproc, dtset%nscforder) 320 ! Shallow copy of the kernel (still owned by ABINIT). 321 wvl%den%denspot%pkernel = kernel 322 wvl%den%denspot%pkernelseq = kernel 323 !Associate the denspot distribution into mpi_enreg. 324 mpi_enreg%nscatterarr => wvl%den%denspot%dpbox%nscatterarr 325 mpi_enreg%ngatherarr => wvl%den%denspot%dpbox%ngatherarr 326 mpi_enreg%ngfft3_ionic = wvl%den%denspot%dpbox%n3pi 327 call wvl_setngfft(me, dtset%mgfft, dtset%nfft, & 328 & dtset%ngfft, nproc, wvl%den%denspot%dpbox%ndims(1), & 329 & wvl%den%denspot%dpbox%ndims(2), & 330 & wvl%den%denspot%dpbox%ndims(3),wvl%den%denspot%dpbox%n3d) 331 332 !We copy the geometry structure. 333 call wvl_wfs_lr_copy(wvl%wfs, wvl%descr) 334 !Reallocate them with new size. 335 call wvl_wfs_set(dtset%strprecon,dtset%spinmagntarget, dtset%kpt, me, dtset%natom, sum(dtset%nband), dtset%nkpt, & 336 & nproc, dtset%nspinor, dtset%nsppol, dtset%nwfshist, dtset%occ_orig(:,1), psps, rprimd, & 337 & wvl%wfs, dtset%wtk, wvl%descr, dtset%wvl_crmult, dtset%wvl_frmult, xred) 338 339 !Recopy old eval for precond. 340 wvl%wfs%ks%orbs%eval = eigen_old 341 ABI_FREE(eigen_old) 342 343 !We allocate psi. 344 !ABI_MALLOC(wvl%wfs%ks%psi,( max(wvl%wfs%ks%orbs%npsidim_comp,wvl%wfs%ks%orbs%npsidim_orbs)+ndebug) ) 345 wvl%wfs%ks%psi=f_malloc_ptr(max(wvl%wfs%ks%orbs%npsidim_comp,wvl%wfs%ks%orbs%npsidim_orbs)+ndebug,id='psi') 346 write(message, '(a,a,a,a,I0)' ) ch10, & 347 & ' wvl_wfsinp_reformat: allocate wavefunctions,', ch10, & 348 & ' size of the compressed array per proc: ', & 349 & product(shape(wvl%wfs%ks%psi)) 350 call wrtout(std_out,message,'COLL') 351 352 !Convert input xred (reduced coordinates) to xcart (cartesian) 353 ABI_MALLOC(xcart,(3, dtset%natom)) 354 call xred2xcart(dtset%natom, rprimd, xcart, xred) 355 356 !We transfer the old wavefunctions to the new ones. 357 call reformatmywaves(me, wvl%wfs%ks%orbs, wvl%descr%atoms, & 358 & hgrid_old(1), hgrid_old(2), hgrid_old(3), nSize_old(1), nSize_old(2), & 359 & nSize_old(3), xcart_old, keys_old, psi_old, wvl%descr%h(1), wvl%descr%h(2), & 360 & wvl%descr%h(3), wvl%descr%Glr%d%n1, wvl%descr%Glr%d%n2, wvl%descr%Glr%d%n3, xcart, & 361 & wvl%wfs%ks%lzd%Glr%wfd, wvl%wfs%ks%psi) 362 ABI_FREE(xcart) 363 ABI_FREE(xcart_old) 364 365 !We free the old descriptors and arrays. 366 ABI_FREE(psi_old) 367 call deallocate_wfd(keys_old) 368 369 call local_potential_dimensions(me,wvl%wfs%ks%lzd,wvl%wfs%ks%orbs,wvl%den%denspot%xc,& 370 & wvl%den%denspot%dpbox%ngatherarr(0,1)) 371 372 !it seems that the table "wvl%projectors%G" is no more used 373 !but it's not allocated -> fortran runtime error 374 #if defined HAVE_BIGDFT 375 ABI_MALLOC(wvl%projectors%G,(dtset%ntypat)) 376 do itypat=1,dtset%ntypat 377 call nullify_gaussian_basis(wvl%projectors%G(itypat)) 378 end do 379 #endif 380 381 !Reallocate projectors for the new positions. 382 call wvl_projectors_set(me, dtset%natom, wvl%projectors, psps, rprimd, & 383 & wvl%wfs, wvl%descr, dtset%wvl_frmult, xred) 384 385 !Orthogonilise new wavefunctions. 386 call first_orthon(me, nproc, wvl%wfs%ks%orbs, wvl%wfs%ks%lzd, wvl%wfs%ks%comms, & 387 & wvl%wfs%ks%psi, wvl%wfs%ks%hpsi, wvl%wfs%ks%psit, wvl%wfs%ks%orthpar,wvl%descr%paw) 388 389 #else 390 BIGDFT_NOTENABLED_ERROR() 391 if (.false.) write(std_out,*) dtset%nstep,mpi_enreg%me,psps%npsp,wvl%wfs%ks,rprimd(1,1),& 392 & xred_old(1,1),xred(1,1) 393 #endif 394 395 end subroutine wvl_wfsinp_reformat
ABINIT/wvl_wfsinp_scratch [ 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=information 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 information on wf file. xred(3,natom)=reduced dimensionless atomic coordinates
OUTPUT
SIDE EFFECTS
wvl <type(wvl_data)>=wavefunctions & projectors information for wavelets.
SOURCE
429 subroutine wvl_wfsinp_scratch(dtset, mpi_enreg, occ, rprimd, wvl, xred) 430 431 #if defined HAVE_BIGDFT 432 use BigDFT_API, only : createIonicPotential, input_wf_diag, gaussian_basis, & 433 & input_variables, calculate_rhocore, deallocate_Lzd_except_Glr, INPUT_IG_OFF,& 434 & SMEARING_DIST_ERF, PSPCODE_PAW 435 #endif 436 437 !Arguments ------------------------------- 438 !scalars 439 type(dataset_type), intent(in) :: dtset 440 type(MPI_type), intent(inout) :: mpi_enreg 441 type(wvl_data), intent(inout) :: wvl 442 !arrays 443 real(dp), intent(inout) :: occ(dtset%mband*dtset%nkpt*dtset%nsppol) 444 real(dp), intent(in) :: rprimd(3, 3) 445 real(dp), intent(in) :: xred(3, dtset%natom) 446 447 !Local variables------------------------------- 448 #if defined HAVE_BIGDFT 449 character(len = 500) :: message 450 integer :: comm,me,nproc 451 integer :: iscf_local 452 integer :: nvirt 453 integer :: ii,shift_vpsp,size_vpsp 454 logical :: onlywf=.false. ! find the wavefunctions and return 455 logical :: wvlbigdft=.false. 456 real(dp), allocatable :: xcart(:,:) 457 real(dp), allocatable :: rhor(:,:) 458 real(dp), pointer :: vpsp(:) 459 real(dp):: elecfield(3) 460 type(gaussian_basis) :: Gvirt 461 type(input_variables) :: in ! To be removed, waiting for BigDFT upgrade 462 #endif 463 464 ! ********************************************************************* 465 466 #if defined HAVE_BIGDFT 467 468 elecfield=zero 469 470 !If usewvl: wvlbigdft indicates that the BigDFT workflow will be followed 471 wvlbigdft=(dtset%usewvl==1.and.dtset%wvl_bigdft_comp==1) 472 473 write(message, '(a,a)' ) ch10,& 474 & ' wvl_wfsinp_scratch: wavefunction initialisation.' 475 call wrtout(std_out,message,'COLL') 476 477 comm=mpi_enreg%comm_wvl 478 me=xmpi_comm_rank(comm) 479 nproc=xmpi_comm_size(comm) 480 !Store xcart for each atom 481 ABI_MALLOC(xcart,(3, dtset%natom)) 482 call xred2xcart(dtset%natom, rprimd, xcart, xred) 483 484 !We allocate temporary arrays for rho and vpsp. 485 !allocate ionic potential 486 if (wvl%den%denspot%dpbox%n3pi > 0) then 487 size_vpsp=wvl%descr%Glr%d%n1i*wvl%descr%Glr%d%n2i*wvl%den%denspot%dpbox%n3pi 488 shift_vpsp=wvl%den%denspot%dpbox%ndims(1)*wvl%den%denspot%dpbox%ndims(2) & 489 & *wvl%den%denspot%dpbox%nscatterarr(me,4) 490 ABI_MALLOC(vpsp,(size_vpsp+shift_vpsp)) 491 else 492 ABI_MALLOC(vpsp,(1)) 493 end if 494 495 if(.not. onlywf) then 496 ABI_MALLOC(rhor,(dtset%nfft,dtset%nspden)) 497 call mklocl_wavelets(elecfield, xcart, mpi_enreg, dtset%natom, & 498 & dtset%nfft, dtset%nspden, 1, rprimd, vpsp, & 499 & wvl%den, wvl%descr, xcart) 500 ABI_FREE(rhor) 501 end if 502 503 ! IMPORTANT: onlywf=.true. does not work yet, do not change this: 504 ! if(.not. wvlbigdft) onlywf=.true. !do not apply Hamiltonian inside input_wf_diag. 505 if(dtset%usepaw==1) wvl%descr%atoms%npspcode(:)=wvl%descr%npspcode_paw_init_guess(:) 506 iscf_local=dtset%iscf 507 if(.not. wvlbigdft) iscf_local=0 !important to have good occ values 508 509 !This routine allocates psi, hpsi and psit. 510 nvirt = 0 511 in%linear = INPUT_IG_OFF 512 in%nspin = dtset%nsppol 513 in%exctxpar = wvl%descr%exctxpar 514 in%itrpmax = dtset%nnsclo 515 in%iscf = iscf_local !dtset%iscf 516 in%Tel = dtset%tsmear 517 ! if (dtset%iscf == 0) in%Tel = zero 518 if (iscf_local == 0) in%Tel = zero 519 in%SIC = wvl%wfs%ks%SIC 520 in%orthpar = wvl%wfs%ks%orthpar 521 522 in%occopt=dtset%occopt 523 if(dtset%occopt>2 .and. dtset%occopt<7)then 524 call wvl_occopt_abi2big(in%occopt,in%occopt,1) 525 else 526 ! This will be used only for the initial wavefunctions: 527 in%occopt=SMEARING_DIST_ERF 528 end if 529 530 !Note: check if all required in "in" is passed. 531 !remove argument "in" from input_wf_diag 532 call input_wf_diag(me, nproc, & 533 & wvl%descr%atoms, wvl%den%denspot,& 534 & wvl%wfs%ks%orbs, nvirt, wvl%wfs%ks%comms, & 535 & wvl%wfs%ks%lzd, wvl%e%energs, xcart, & 536 & wvl%projectors%nlpsp, & 537 & dtset%ixc, wvl%wfs%ks%psi, wvl%wfs%ks%hpsi, wvl%wfs%ks%psit, & 538 & Gvirt, in%nspin, wvl%wfs%GPU, in, onlywf, wvl%projectors%G, wvl%descr%paw) 539 540 !This provisory: wvl%descr%paw could be passed as optional 541 !to input_wf_diag to allocate spsi inside this routine 542 if(dtset%usepaw==1) then 543 wvl%descr%atoms%npspcode(:)=PSPCODE_PAW 544 ABI_MALLOC(wvl%descr%paw%spsi,(max(wvl%wfs%ks%orbs%npsidim_orbs,wvl%wfs%ks%orbs%npsidim_comp))) 545 do ii=1,size(wvl%wfs%ks%psi) 546 wvl%descr%paw%spsi(ii)=wvl%wfs%ks%psi(ii) 547 wvl%wfs%ks%hpsi(ii)=wvl%wfs%ks%psi(ii) 548 end do 549 end if 550 551 if(wvlbigdft ) then 552 ! Copy occupations from BigDFT objects to ABINIT 553 call wvl_occ_abi2big(dtset%mband,dtset%nkpt,dtset%nsppol,occ,2,wvl%wfs) 554 else 555 ! Copy occupations from ABINIT to BigDFT objects 556 call wvl_occ_abi2big(dtset%mband,dtset%nkpt,dtset%nsppol,occ,1,wvl%wfs) 557 end if 558 559 write(message, '(a)' ) & 560 & ' | wavefunctions have been calculated.' 561 call wrtout(std_out,message,'COLL') 562 563 ABI_FREE(xcart) 564 ABI_FREE(vpsp) 565 566 #else 567 BIGDFT_NOTENABLED_ERROR() 568 if (.false.) write(std_out,*) dtset%nstep,mpi_enreg%me,wvl%wfs%ks,occ(1),rprimd(1,1),xred(1,1) 569 #endif 570 571 end subroutine wvl_wfsinp_scratch