TABLE OF CONTENTS
ABINIT/wvl_read [ Functions ]
NAME
wvl_read
FUNCTION
Simple wrapper around the read disk methods of BigDFT for 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 . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
dtset <type(dataset_type)>=input variables. mpi_enreg=informations about MPI parallelization option= -2 for write with BigDFT format, -1 for reading wavelets coefficients with BigDFT format, 2 for write, 1 for read. rprimd(3,3)=dimensional primitive translations in real space (bohr) wff <type(wffile_type)>=struct info for wavefunction wfs <type(wvl_wf_type)>=wavefunctions informations for wavelets. xred(3,natom)=reduced dimensionless atomic coordinates
OUTPUT
PARENTS
wvl_wfsinp_disk
CHILDREN
etsf_io_basisdata_put,etsf_io_electrons_put,etsf_io_main_put writemywaves,writeonewave,wrtout,xred2xcart
SOURCE
39 #if defined HAVE_CONFIG_H 40 #include "config.h" 41 #endif 42 43 #include "abi_common.h" 44 45 subroutine wvl_read(dtset, hdr0, hdr, mpi_enreg, option, rprimd, wff, wfs, wvl, xred) 46 47 use defs_basis 48 use defs_abitypes 49 use defs_abitypes 50 use defs_wvltypes 51 use m_wffile 52 use m_errors 53 use m_profiling_abi 54 use m_xmpi 55 56 use m_geometry, only : xred2xcart 57 #if defined HAVE_BIGDFT 58 use BigDFT_API, only: readonewave, reformatonewave, readmywaves, & 59 & WF_FORMAT_NONE 60 #endif 61 #if defined HAVE_ETSF_IO 62 use etsf_io 63 #endif 64 65 !This section has been created automatically by the script Abilint (TD). 66 !Do not modify the following lines by hand. 67 #undef ABI_FUNC 68 #define ABI_FUNC 'wvl_read' 69 use interfaces_14_hidewrite 70 !End of the abilint section 71 72 implicit none 73 74 !Arguments ------------------------------- 75 !scalars 76 integer, intent(in) :: option 77 type(dataset_type), intent(in) :: dtset 78 type(hdr_type), intent(in) :: hdr0 79 type(hdr_type), intent(in) :: hdr 80 type(MPI_type), intent(in) :: mpi_enreg 81 type(wffile_type),intent(in) :: wff 82 type(wvl_wf_type), intent(inout) :: wfs 83 type(wvl_internal_type), intent(in) :: wvl 84 !arrays 85 real(dp), intent(in) :: rprimd(3, 3) 86 real(dp), intent(in) :: xred(3, dtset%natom) 87 88 !Local variables------------------------------- 89 #if defined HAVE_BIGDFT 90 character(len = 500) :: message 91 integer :: iBand, bandSize 92 integer :: comm,me 93 real(dp), allocatable :: xcart(:,:), psifscf(:,:,:) 94 real(dp), allocatable :: xcart_old(:,:) 95 #if defined HAVE_ETSF_IO 96 integer :: i, iCoeff, iCoarse, iFine 97 integer :: n_old(3) 98 real(dp) :: hgrid_old(3) 99 real(dp), allocatable :: psigold(:,:,:,:,:,:) 100 real(dp), allocatable, target :: wvl_band(:) 101 type(etsf_basisdata) :: basis_folder 102 type(etsf_main) :: main_folder 103 type(etsf_electrons) :: electrons_folder 104 logical :: reformat 105 logical :: lstat 106 type(etsf_io_low_error) :: error 107 #endif 108 #endif 109 110 ! ********************************************************************* 111 112 #if defined HAVE_BIGDFT 113 114 if (abs(option) /= 1) then 115 write(message,'(a,a,a,i0,a)')& 116 & ' Option argument is wrong,', ch10, & 117 & ' awaited values are -1 or 1 but option = ', option, '.' 118 MSG_BUG(message) 119 end if 120 121 comm=mpi_enreg%comm_wvl 122 me=xmpi_comm_rank(comm) 123 !Store xcart for each atom 124 ABI_ALLOCATE(xcart,(3, dtset%natom)) 125 ABI_ALLOCATE(xcart_old,(3, dtset%natom)) 126 call xred2xcart(dtset%natom, rprimd, xcart, xred) 127 128 write(message,'(2a)') ch10,' wvl_read: read wavefunctions from file.' 129 call wrtout(std_out,message,'COLL') 130 131 if (option > 0) then 132 bandSize = wfs%ks%lzd%Glr%wfd%nvctr_c + 7 * wfs%ks%lzd%Glr%wfd%nvctr_f 133 ! Read in the ABINIT way. 134 if (wff%iomode == IO_MODE_FORTRAN .or. (wff%iomode == IO_MODE_FORTRAN_MASTER .and. wff%master==wff%me)) then 135 ABI_ALLOCATE(psifscf,(wvl%Glr%d%n1i, wvl%Glr%d%n2i, wvl%Glr%d%n3i)) 136 do iBand = 1, dtset%mband * dtset%nsppol, 1 137 call readonewave(wff%unwff, .false., iBand, me, & 138 & wvl%Glr%d%n1, wvl%Glr%d%n2, wvl%Glr%d%n3, & 139 & wvl%h(1), wvl%h(2), wvl%h(3), wvl%atoms, & 140 & wfs%ks%lzd%Glr%wfd, xcart_old, xcart, & 141 & wfs%ks%psi(bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + 1: & 142 & bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + bandSize), & 143 & wfs%ks%orbs%eval(iBand), psifscf) 144 end do 145 ABI_DEALLOCATE(psifscf) 146 #if defined HAVE_ETSF_IO 147 else if (wff%iomode == IO_MODE_ETSF) then 148 ! Read a NetCDF file. 149 ! coordinates_of_grid_points is used to store the geometric 150 ! position of coefficients of wavelets i, as integer in 151 ! wvl%ni(:) dimensions. 152 ! coefficients_of_wavefunctions is used to store the psi values for 153 ! each wavelet. 154 155 ! We check if we need reformating 156 if (abs(hdr%rprimd(1,1) / hdr%ngfft(1) - & 157 & hdr0%rprimd(1,1) / hdr0%ngfft(1)) < tol8 .and. & 158 & maxval(abs(hdr%nwvlarr - hdr0%nwvlarr)) == 0 .and. & 159 & maxval(abs(hdr%ngfft - hdr0%ngfft )) == 0 ) then 160 reformat = .false. 161 write(message, '(2a)' ) ch10,& 162 & ' wvl_read: wavefunctions needs NO reformatting.' 163 call wrtout(std_out,message,'COLL') 164 165 else 166 reformat = .true. 167 write(message, '(2a)' ) ch10,& 168 & ' wvl_read: wavefunctions needs reformatting.' 169 call wrtout(std_out,message,'COLL') 170 171 ! Values for old system 172 call xred2xcart(dtset%natom, hdr0%rprimd, xcart_old, hdr0%xred) 173 displ = zero 174 do i = 1, dtset%natom, 1 175 displ = (xcart_old(1, i) - xcart(1, i)) ** 2 + & 176 & (xcart_old(2, i) - xcart(2, i)) ** 2 + & 177 & (xcart_old(3, i) - xcart(3, i)) ** 2 178 end do 179 displ = sqrt(displ) 180 181 ! Do the reformating by expanding old wavefunctions on the grid 182 ! and interpolate it on the new grid using BigDFT reformatonewave(). 183 n_old = (hdr0%ngfft(1:3) - 31) / 2 184 hgrid_old(1) = hdr0%rprimd(1,1) / n_old(1) 185 hgrid_old(2) = hdr0%rprimd(2,2) / n_old(2) 186 hgrid_old(3) = hdr0%rprimd(3,3) / n_old(3) 187 ABI_ALLOCATE(psigold,(0:n_old(1), 2, 0:n_old(2), 2, 0:n_old(3), 2)) 188 ABI_ALLOCATE(psifscf,(wvl%Glr%d%n1i, wvl%Glr%d%n2i, wvl%Glr%d%n3i)) 189 ABI_ALLOCATE(basis_folder%coordinates_of_basis_grid_points%data2D,(3,hdr0%nwvlarr(1))) 190 end if 191 192 ! We read the basis set definition 193 ABI_ALLOCATE(basis_folder%number_of_coefficients_per_grid_point%data1D,(hdr0%nwvlarr(1))) 194 call etsf_io_basisdata_get(wff%unwff, basis_folder, lstat, error) 195 ETSF_CHECK_ERROR(lstat, error) 196 197 ! We allocate a temporary array to read the wavefunctions 198 ! and reorder then into wfs%ks%psi 199 ABI_ALLOCATE(wvl_band,(sum(basis_folder%number_of_coefficients_per_grid_point%data1D))) 200 ! We just read the file from disk to memory, band per band. 201 do iBand = 1, dtset%mband, 1 202 main_folder%wfs_coeff__state_access = iBand 203 main_folder%coefficients_of_wavefunctions%data1D => wvl_band 204 call etsf_io_main_get(wff%unwff, main_folder, lstat, error) 205 ETSF_CHECK_ERROR(lstat, error) 206 207 ! Now we reorder 208 if (reformat) then 209 ! Expand the wf on the old grid. 210 psigold = zero 211 iCoeff = 1 212 do i = 1, hdr0%nwvlarr(1), 1 213 coord = basis_folder%coordinates_of_basis_grid_points%data2D(:, i) / 2 214 psigold(coord(1), 1, coord(2), 1, coord(3), 1) = wvl_band(iCoeff) 215 iCoeff = iCoeff + 1 216 if (basis_folder%number_of_coefficients_per_grid_point%data1D(i) == 8) then 217 psigold(coord(1), 2, coord(2), 1, coord(3), 1) = wvl_band(iCoeff + 0) 218 psigold(coord(1), 1, coord(2), 2, coord(3), 1) = wvl_band(iCoeff + 1) 219 psigold(coord(1), 2, coord(2), 2, coord(3), 1) = wvl_band(iCoeff + 2) 220 psigold(coord(1), 1, coord(2), 1, coord(3), 2) = wvl_band(iCoeff + 3) 221 psigold(coord(1), 2, coord(2), 1, coord(3), 2) = wvl_band(iCoeff + 4) 222 psigold(coord(1), 1, coord(2), 2, coord(3), 2) = wvl_band(iCoeff + 5) 223 psigold(coord(1), 2, coord(2), 2, coord(3), 2) = wvl_band(iCoeff + 6) 224 iCoeff = iCoeff + 7 225 end if 226 end do 227 228 call reformatonewave(displ, wfs%ks%lzd%Glr%wfd, wvl%atoms, hgrid_old(1), & 229 & hgrid_old(2), hgrid_old(3), n_old(1), n_old(2), & 230 & n_old(3), xcart_old, psigold, & 231 & wvl%h(1), wvl%h(2), wvl%h(3), & 232 & wvl%Glr%d%n1, wvl%Glr%d%n2, & 233 & wvl%Glr%d%n3, xcart, psifscf, & 234 & wfs%ks%psi(bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + 1: & 235 & bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + bandSize)) 236 else 237 ! No reformating, just reorder 238 iCoeff = 1 239 iFine = bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + wfs%ks%lzd%Glr%wfd%nvctr_c + 1 240 do iCoarse = 1, hdr0%nwvlarr(1), 1 241 wfs%ks%psi(bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + iCoarse) = & 242 & wvl_band(iCoeff) 243 iCoeff = iCoeff + 1 244 if (basis_folder%number_of_coefficients_per_grid_point%data1D(iCoarse) == 8) then 245 wfs%ks%psi(iFine:iFine + 6) = wvl_band(iCoeff:iCoeff + 6) 246 iFine = iFine + 7 247 iCoeff = iCoeff + 7 248 end if 249 end do 250 end if 251 end do 252 253 ! We read wfs%ks%eval (TODO maybe removed later). 254 electrons_folder%eigenvalues%data1D => wfs%ks%orbs%eval 255 call etsf_io_electrons_get(wff%unwff, electrons_folder, lstat, error) 256 ETSF_CHECK_ERROR(lstat, error) 257 258 ! Deallocate temporary arrays 259 ABI_DEALLOCATE(wvl_band) 260 ABI_DEALLOCATE(basis_folder%number_of_coefficients_per_grid_point%data1D) 261 if (reformat) then 262 ABI_DEALLOCATE(basis_folder%coordinates_of_basis_grid_points%data2D) 263 ABI_DEALLOCATE(psigold) 264 ABI_DEALLOCATE(psifscf) 265 end if 266 #endif 267 268 else 269 write(message,'(4a,i0,a)') ch10,& 270 & ' wff%iomode argument is wrong,', ch10, & 271 & ' awaited values are -1, 0 (or 3 if netcdf/etsf_io is available) but value = ', wff%iomode, '.' 272 MSG_BUG(message) 273 end if 274 else 275 call readmywaves(me, "wavefunctions", WF_FORMAT_NONE, wfs%ks%orbs, & 276 & wvl%Glr%d%n1, wvl%Glr%d%n2, wvl%Glr%d%n3, & 277 & wvl%h(1), wvl%h(2), wvl%h(3), wvl%atoms, & 278 & xcart_old, xcart, wfs%ks%lzd%Glr%wfd, wfs%ks%psi) 279 end if 280 281 ABI_DEALLOCATE(xcart) 282 ABI_DEALLOCATE(xcart_old) 283 #else 284 BIGDFT_NOTENABLED_ERROR() 285 if (.false.) write(std_out,*) option,dtset%nstep,hdr0%ecut,hdr%ecut,mpi_enreg%nproc,wff%me,& 286 & wfs%ks,wvl%h(1),rprimd(1,1),xred(1,1) 287 #endif 288 289 end subroutine wvl_read
ABINIT/wvl_write [ Functions ]
NAME
wvl_write
FUNCTION
Simple wrapper around the write disk methods of BigDFT for wavefunctions.
INPUTS
dtset <type(dataset_type)>=input variables. mpi_enreg=informations about MPI parallelization option= -2 for write with BigDFT format, -1 for reading wavelets coefficients with BigDFT format, 2 for write, 1 for read. rprimd(3,3)=dimensional primitive translations in real space (bohr) wff <type(wffile_type)>=struct info for wavefunction wfs <type(wvl_wf_type)>=wavefunctions informations for wavelets.
OUTPUT
SIDE EFFECTS
xred(3,natom)=reduced dimensionless atomic coordinates
PARENTS
m_iowf
CHILDREN
etsf_io_basisdata_put,etsf_io_electrons_put,etsf_io_main_put writemywaves,writeonewave,wrtout,xred2xcart
SOURCE
324 subroutine wvl_write(dtset, eigen, mpi_enreg, option, rprimd, wff, wfs, wvl, xred) 325 326 use defs_basis 327 use defs_abitypes 328 use defs_wvltypes 329 use m_wffile 330 use m_errors 331 use m_profiling_abi 332 use m_xmpi 333 use m_geometry, only : xred2xcart 334 335 #if defined HAVE_BIGDFT 336 use BigDFT_API, only : writeonewave,writemywaves,WF_FORMAT_NONE 337 #endif 338 #if defined HAVE_ETSF_IO 339 use etsf_io 340 #endif 341 342 !This section has been created automatically by the script Abilint (TD). 343 !Do not modify the following lines by hand. 344 #undef ABI_FUNC 345 #define ABI_FUNC 'wvl_write' 346 use interfaces_14_hidewrite 347 !End of the abilint section 348 349 implicit none 350 351 !Arguments ------------------------------- 352 !scalars 353 integer, intent(in) :: option 354 type(dataset_type), intent(in) :: dtset 355 type(MPI_type), intent(in) :: mpi_enreg 356 type(wffile_type),intent(in) :: wff 357 type(wvl_wf_type), intent(in) :: wfs 358 type(wvl_internal_type), intent(in) :: wvl 359 !arrays 360 real(dp), intent(in), target :: eigen(dtset%mband) 361 real(dp), intent(in) :: rprimd(3, 3) 362 real(dp), intent(in) :: xred(3, dtset%natom) 363 364 !Local variables------------------------------- 365 #if defined HAVE_BIGDFT 366 character(len = 500) :: message 367 integer :: comm,me 368 integer :: iorb 369 integer :: iseg, nseg, ipsi, npsi 370 real(dp), allocatable :: xcart(:,:) 371 #if defined HAVE_ETSF_IO 372 integer :: i, i0, i1, i2, i3, jj, j0, j1, ii 373 integer :: iGrid, iCoeff, iCoarse, iFine 374 integer :: iband 375 integer, allocatable :: coeff_map(:,:,:) 376 integer, allocatable, target :: wvl_coord(:,:), wvl_ncoeff(:) 377 real(dp), allocatable, target :: wvl_band(:) 378 type(etsf_basisdata) :: basis_folder 379 type(etsf_main) :: main_folder 380 type(etsf_electrons) :: electrons_folder 381 logical :: lstat 382 type(etsf_io_low_error) :: error 383 #endif 384 #endif 385 386 ! ********************************************************************* 387 388 #if defined HAVE_BIGDFT 389 390 if (abs(option) /= 2) then 391 write(message,'(a,a,a,i0,a)')& 392 & ' Option argument is wrong,', ch10, & 393 & ' awaited values are -2 or 2 but option = ', option, '.' 394 MSG_BUG(message) 395 end if 396 397 comm=mpi_enreg%comm_wvl 398 me=xmpi_comm_rank(comm) 399 !Store xcart for each atom 400 ABI_ALLOCATE(xcart,(3, dtset%natom)) 401 call xred2xcart(dtset%natom, rprimd, xcart, xred) 402 403 write(message, '(a,a,a,a)' ) ch10,& 404 & ' wvl_write: Write wavefunctions to file.' 405 call wrtout(std_out,message,'COLL') 406 407 if (option > 0) then 408 ! Write in the ABINIT way. 409 if (wff%iomode == IO_MODE_FORTRAN .or. (wff%iomode == IO_MODE_FORTRAN_MASTER .and. wff%master==wff%me)) then 410 iseg = wfs%ks%lzd%Glr%wfd%nseg_c 411 nseg = wfs%ks%lzd%Glr%wfd%nseg_c + wfs%ks%lzd%Glr%wfd%nseg_f 412 ipsi = wfs%ks%lzd%Glr%wfd%nvctr_c 413 npsi = wfs%ks%lzd%Glr%wfd%nvctr_c + 7 * wfs%ks%lzd%Glr%wfd%nvctr_f 414 do iorb = 1, dtset%mband 415 call writeonewave(wff%unwff, .false., iorb, wvl%Glr%d%n1, & 416 & wvl%Glr%d%n2, wvl%Glr%d%n3, & 417 & wvl%h(1), wvl%h(2), wvl%h(3), dtset%natom, & 418 & xcart, wfs%ks%lzd%Glr%wfd%nseg_c, wfs%ks%lzd%Glr%wfd%nvctr_c, & 419 & wfs%ks%lzd%Glr%wfd%keygloc(:,1:iseg), wfs%ks%lzd%Glr%wfd%keyvloc(1:iseg), wfs%ks%lzd%Glr%wfd%nseg_f, & 420 & wfs%ks%lzd%Glr%wfd%nvctr_f, wfs%ks%lzd%Glr%wfd%keygloc(:, iseg + 1:nseg), & 421 & wfs%ks%lzd%Glr%wfd%keyvloc(iseg + 1:nseg), & 422 & wfs%ks%psi(npsi * (iorb - me * wfs%ks%orbs%norbp - 1) + 1: & 423 & npsi * (iorb - me * wfs%ks%orbs%norbp - 1) + ipsi), & 424 & wfs%ks%psi(npsi * (iorb - me * wfs%ks%orbs%norbp - 1) + ipsi + 1: & 425 & npsi * (iorb - me * wfs%ks%orbs%norbp - 1) + npsi), & 426 & wfs%ks%orbs%eval(iorb)) 427 end do 428 #if defined HAVE_ETSF_IO 429 else if (wff%iomode == IO_MODE_ETSF) then 430 ! Write a NetCDF file. 431 ! coordinates_of_grid_points is used to store the geometric 432 ! position of coefficients of wavelets i, as integer in 433 ! wvl%dpSize(:) dimensions. 434 ! coefficients_of_wavefunctions is used to store the psi values for 435 ! each wavelet. 436 437 ! We write the basis set. 438 ! ======================= 439 ABI_ALLOCATE(wvl_coord,(3, wfs%ks%lzd%Glr%wfd%nvctr_c)) 440 ABI_ALLOCATE(wvl_ncoeff,(wfs%ks%lzd%Glr%wfd%nvctr_c)) 441 ! Will store the grid index for a given geometric point 442 ABI_ALLOCATE(coeff_map,(0:wvl%Glr%d%n1,0:wvl%Glr%d%n2,0:wvl%Glr%d%n3)) 443 ! coarse part 444 iGrid = 0 445 coeff_map = 0 446 do iseg = 1, wfs%ks%lzd%Glr%wfd%nseg_c 447 jj = wfs%ks%lzd%Glr%wfd%keyvloc(iseg) 448 j0 = wfs%ks%lzd%Glr%wfd%keygloc(1, iseg) 449 j1 = wfs%ks%lzd%Glr%wfd%keygloc(2, iseg) 450 ii = j0 - 1 451 i3 = ii / ((wvl%Glr%d%n1 + 1) * & 452 & (wvl%Glr%d%n2 + 1)) 453 ii = ii - i3 * (wvl%Glr%d%n1 + 1) * & 454 & (wvl%Glr%d%n2 + 1) 455 i2 = ii / (wvl%Glr%d%n1 + 1) 456 i0 = ii - i2 * (wvl%Glr%d%n1 + 1) 457 i1 = i0 + j1 - j0 458 do i = i0, i1 459 iGrid = iGrid + 1 460 coeff_map(i, i2, i3) = iGrid 461 wvl_coord(:, iGrid) = (/ i * 2, i2 * 2, i3 * 2 /) 462 wvl_ncoeff(iGrid) = 1 463 end do 464 end do 465 ! fine part 466 do iseg = 1, wfs%ks%lzd%Glr%wfd%nseg_f 467 jj = wfs%ks%lzd%Glr%wfd%keyvloc(wfs%ks%lzd%Glr%wfd%nseg_c + iseg) 468 j0 = wfs%ks%lzd%Glr%wfd%keygloc(1, wfs%ks%lzd%Glr%wfd%nseg_c + iseg) 469 j1 = wfs%ks%lzd%Glr%wfd%keygloc(2, wfs%ks%lzd%Glr%wfd%nseg_c + iseg) 470 ii = j0 - 1 471 i3 = ii / ((wvl%Glr%d%n1 + 1) * & 472 & (wvl%Glr%d%n2 + 1)) 473 ii = ii - i3 * (wvl%Glr%d%n1 + 1) * & 474 & (wvl%Glr%d%n2 + 1) 475 i2 = ii / (wvl%Glr%d%n1 + 1) 476 i0 = ii - i2 * (wvl%Glr%d%n1 + 1) 477 i1 = i0 + j1 - j0 478 do i = i0, i1 479 iGrid = coeff_map(i , i2 , i3) 480 wvl_ncoeff(iGrid) = wvl_ncoeff(iGrid) + 7 481 end do 482 end do 483 basis_folder%coordinates_of_basis_grid_points%data2D => wvl_coord 484 basis_folder%number_of_coefficients_per_grid_point%data1D => wvl_ncoeff 485 486 call etsf_io_basisdata_put(wff%unwff, basis_folder, lstat, error) 487 488 ETSF_CHECK_ERROR(lstat, error) 489 490 ! We write the wavefuctions. 491 ! ========================== 492 ! We write the wavefunctions band per band. 493 ABI_ALLOCATE(wvl_band,(wfs%ks%lzd%Glr%wfd%nvctr_c + 7 * wfs%ks%lzd%Glr%wfd%nvctr_f)) 494 ! We reorder the wfs%ks%psi(:,iband) into wvl_band 495 do iband = 1, dtset%mband, 1 496 ! iCoeff is the index of the coefficient we are writing in wvl_band 497 iCoeff = 1 498 ! iCoarse is the index of the coarse part in wfs%ks%psi (=iGrid) 499 iCoarse = 1 500 ! iFine is the index of the fine part in wfs%ks%psi 501 iFine = wfs%ks%lzd%Glr%wfd%nvctr_c + 1 502 npsi = wfs%ks%lzd%Glr%wfd%nvctr_c + 7 * wfs%ks%lzd%Glr%wfd%nvctr_f 503 do iGrid = 1, wfs%ks%lzd%Glr%wfd%nvctr_c, 1 504 wvl_band(iCoeff) = & 505 & wfs%ks%psi((iband - me * wfs%ks%orbs%norbp - 1) * npsi + iCoarse) 506 iCoeff = iCoeff + 1 507 iCoarse = iCoarse + 1 508 if (wvl_ncoeff(iGrid) == 8) then 509 wvl_band(iCoeff:iCoeff + 6) = & 510 & wfs%ks%psi((iband - me * wfs%ks%orbs%norbp - 1) * npsi + iFine: & 511 & (iband - me * wfs%ks%orbs%norbp - 1) * npsi + iFine + 6) 512 iCoeff = iCoeff + 7 513 iFine = iFine + 7 514 end if 515 end do 516 main_folder%wfs_coeff__state_access = iband 517 main_folder%coefficients_of_wavefunctions%data1D => wvl_band 518 call etsf_io_main_put(wff%unwff, main_folder, lstat, error) 519 ETSF_CHECK_ERROR(lstat, error) 520 end do 521 522 ! We write the electronic informations. 523 ! ===================================== 524 electrons_folder%eigenvalues%data1D => eigen 525 electrons_folder%eigenvalues__number_of_states = dtset%mband 526 call etsf_io_electrons_put(wff%unwff, electrons_folder, lstat, error) 527 ETSF_CHECK_ERROR(lstat,error) 528 529 ! We deallocate all arrays 530 ABI_DEALLOCATE(wvl_band) 531 ABI_DEALLOCATE(wvl_coord) 532 ABI_DEALLOCATE(wvl_ncoeff) 533 ABI_DEALLOCATE(coeff_map) 534 #endif 535 else 536 write(message,'(3a,i0,a)')& 537 & ' wff%iomode argument is wrong,', ch10, & 538 & ' awaited values are -1, 0 (or 3 if netcdf/etsf_io is available) but value = ', wff%iomode, '.' 539 MSG_BUG(message) 540 end if 541 else 542 ! Write in the BigDFT way. 543 call writemywaves(me, "wavefunctions", WF_FORMAT_NONE, wfs%ks%orbs, & 544 & wvl%Glr%d%n1, wvl%Glr%d%n2, wvl%Glr%d%n3, & 545 & wvl%h(1), wvl%h(2), wvl%h(3),wvl%atoms, & 546 & xcart, wfs%ks%lzd%Glr%wfd, wfs%ks%psi) 547 end if 548 549 ABI_DEALLOCATE(xcart) 550 551 #else 552 BIGDFT_NOTENABLED_ERROR() 553 if (.false.) write(std_out,*) option,dtset%nstep,mpi_enreg%nproc,wff%me,& 554 & wfs%ks,wvl%h(1),eigen(1),rprimd(1,1),xred(1,1) 555 #endif 556 557 end subroutine wvl_write