TABLE OF CONTENTS


ABINIT/m_wvl_rwwf [ Modules ]

[ Top ] [ Modules ]

NAME

  m_wvl_rwwf

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

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_wvl_rwwf
28 
29  use defs_basis
30  use defs_abitypes
31  use defs_datatypes
32  use defs_wvltypes
33  use m_wffile
34  use m_errors
35  use m_abicore
36  use m_xmpi
37 #if defined HAVE_ETSF_IO
38   use etsf_io
39 #endif
40 
41  use m_geometry,     only : xred2xcart
42 
43  implicit none
44 
45  private

ABINIT/wvl_read [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_read

FUNCTION

 Simple wrapper around the read 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.
  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

 85 subroutine wvl_read(dtset, hdr0, hdr, mpi_enreg, option, rprimd, wff, wfs, wvl, xred)
 86 
 87 #if defined HAVE_BIGDFT
 88   use BigDFT_API, only: readonewave, reformatonewave, readmywaves, &
 89 &                       WF_FORMAT_NONE
 90 #endif
 91 
 92 !This section has been created automatically by the script Abilint (TD).
 93 !Do not modify the following lines by hand.
 94 #undef ABI_FUNC
 95 #define ABI_FUNC 'wvl_read'
 96 !End of the abilint section
 97 
 98   implicit none
 99 
100 !Arguments -------------------------------
101   !scalars
102   integer, intent(in)                       :: option
103   type(dataset_type), intent(in)            :: dtset
104   type(hdr_type), intent(in)                :: hdr0
105   type(hdr_type), intent(in)                :: hdr
106   type(MPI_type), intent(in)                :: mpi_enreg
107   type(wffile_type),intent(in)              :: wff
108   type(wvl_wf_type), intent(inout)          :: wfs
109   type(wvl_internal_type), intent(in)       :: wvl
110   !arrays
111   real(dp), intent(in)                      :: rprimd(3, 3)
112   real(dp), intent(in)                      :: xred(3, dtset%natom)
113 
114 !Local variables-------------------------------
115 #if defined HAVE_BIGDFT
116   character(len = 500)  :: message
117   integer               :: iBand, bandSize
118   integer               :: comm,me
119   real(dp), allocatable :: xcart(:,:), psifscf(:,:,:)
120   real(dp), allocatable :: xcart_old(:,:)
121 #if defined HAVE_ETSF_IO
122   integer               :: i, iCoeff, iCoarse, iFine
123   integer               :: n_old(3)
124   real(dp)              :: hgrid_old(3)
125   real(dp), allocatable :: psigold(:,:,:,:,:,:)
126   real(dp), allocatable, target :: wvl_band(:)
127   type(etsf_basisdata) :: basis_folder
128   type(etsf_main)      :: main_folder
129   type(etsf_electrons) :: electrons_folder
130   logical               :: reformat
131   logical              :: lstat
132   type(etsf_io_low_error) :: error
133 #endif
134 #endif
135 
136 ! *********************************************************************
137 
138 #if defined HAVE_BIGDFT
139 
140  if (abs(option) /= 1) then
141    write(message,'(a,a,a,i0,a)')&
142 &   '  Option argument is wrong,', ch10, &
143 &   '  awaited values are -1 or  1 but option = ', option, '.'
144    MSG_BUG(message)
145  end if
146 
147  comm=mpi_enreg%comm_wvl
148  me=xmpi_comm_rank(comm)
149 !Store xcart for each atom
150  ABI_ALLOCATE(xcart,(3, dtset%natom))
151  ABI_ALLOCATE(xcart_old,(3, dtset%natom))
152  call xred2xcart(dtset%natom, rprimd, xcart, xred)
153 
154  write(message,'(2a)') ch10,' wvl_read:  read wavefunctions from file.'
155  call wrtout(std_out,message,'COLL')
156 
157  if (option > 0) then
158    bandSize = wfs%ks%lzd%Glr%wfd%nvctr_c + 7 * wfs%ks%lzd%Glr%wfd%nvctr_f
159 !  Read in the ABINIT way.
160    if (wff%iomode == IO_MODE_FORTRAN .or. (wff%iomode == IO_MODE_FORTRAN_MASTER .and. wff%master==wff%me)) then
161      ABI_ALLOCATE(psifscf,(wvl%Glr%d%n1i, wvl%Glr%d%n2i, wvl%Glr%d%n3i))
162      do iBand = 1, dtset%mband * dtset%nsppol, 1
163        call readonewave(wff%unwff, .false., iBand, me, &
164 &       wvl%Glr%d%n1, wvl%Glr%d%n2, wvl%Glr%d%n3, &
165 &       wvl%h(1), wvl%h(2), wvl%h(3), wvl%atoms, &
166 &       wfs%ks%lzd%Glr%wfd, xcart_old, xcart, &
167 &       wfs%ks%psi(bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + 1: &
168 &       bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + bandSize), &
169 &       wfs%ks%orbs%eval(iBand), psifscf)
170      end do
171      ABI_DEALLOCATE(psifscf)
172 #if defined HAVE_ETSF_IO
173    else if (wff%iomode == IO_MODE_ETSF) then
174 !    Read a NetCDF file.
175 !    coordinates_of_grid_points is used to store the geometric
176 !    position of coefficients of wavelets i, as integer in
177 !    wvl%ni(:) dimensions.
178 !    coefficients_of_wavefunctions is used to store the psi values for
179 !    each wavelet.
180 
181 !    We check if we need reformating
182      if (abs(hdr%rprimd(1,1) / hdr%ngfft(1) - &
183 &     hdr0%rprimd(1,1) / hdr0%ngfft(1)) < tol8 .and. &
184 &     maxval(abs(hdr%nwvlarr - hdr0%nwvlarr)) == 0 .and. &
185 &     maxval(abs(hdr%ngfft   - hdr0%ngfft  )) == 0 ) then
186        reformat = .false.
187        write(message, '(2a)' ) ch10,&
188 &       ' wvl_read:  wavefunctions needs NO reformatting.'
189        call wrtout(std_out,message,'COLL')
190 
191      else
192        reformat = .true.
193        write(message, '(2a)' ) ch10,&
194 &       ' wvl_read:  wavefunctions needs reformatting.'
195        call wrtout(std_out,message,'COLL')
196 
197 !      Values for old system
198        call xred2xcart(dtset%natom, hdr0%rprimd, xcart_old, hdr0%xred)
199        displ = zero
200        do i = 1, dtset%natom, 1
201          displ = (xcart_old(1, i) - xcart(1, i)) ** 2 + &
202 &         (xcart_old(2, i) - xcart(2, i)) ** 2 + &
203 &         (xcart_old(3, i) - xcart(3, i)) ** 2
204        end do
205        displ = sqrt(displ)
206 
207 !      Do the reformating by expanding old wavefunctions on the grid
208 !      and interpolate it on the new grid using BigDFT reformatonewave().
209        n_old = (hdr0%ngfft(1:3) - 31) / 2
210        hgrid_old(1) = hdr0%rprimd(1,1) / n_old(1)
211        hgrid_old(2) = hdr0%rprimd(2,2) / n_old(2)
212        hgrid_old(3) = hdr0%rprimd(3,3) / n_old(3)
213        ABI_ALLOCATE(psigold,(0:n_old(1), 2, 0:n_old(2), 2, 0:n_old(3), 2))
214        ABI_ALLOCATE(psifscf,(wvl%Glr%d%n1i, wvl%Glr%d%n2i, wvl%Glr%d%n3i))
215        ABI_ALLOCATE(basis_folder%coordinates_of_basis_grid_points%data2D,(3,hdr0%nwvlarr(1)))
216      end if
217 
218 !    We read the basis set definition
219      ABI_ALLOCATE(basis_folder%number_of_coefficients_per_grid_point%data1D,(hdr0%nwvlarr(1)))
220      call etsf_io_basisdata_get(wff%unwff, basis_folder, lstat, error)
221      ETSF_CHECK_ERROR(lstat, error)
222 
223 !    We allocate a temporary array to read the wavefunctions
224 !    and reorder then into wfs%ks%psi
225      ABI_ALLOCATE(wvl_band,(sum(basis_folder%number_of_coefficients_per_grid_point%data1D)))
226 !    We just read the file from disk to memory, band per band.
227      do iBand = 1, dtset%mband, 1
228        main_folder%wfs_coeff__state_access = iBand
229        main_folder%coefficients_of_wavefunctions%data1D => wvl_band
230        call etsf_io_main_get(wff%unwff, main_folder, lstat, error)
231        ETSF_CHECK_ERROR(lstat, error)
232 
233 !      Now we reorder
234        if (reformat) then
235 !        Expand the wf on the old grid.
236          psigold = zero
237          iCoeff = 1
238          do i = 1, hdr0%nwvlarr(1), 1
239            coord = basis_folder%coordinates_of_basis_grid_points%data2D(:, i) / 2
240            psigold(coord(1), 1, coord(2), 1, coord(3), 1) = wvl_band(iCoeff)
241            iCoeff = iCoeff + 1
242            if (basis_folder%number_of_coefficients_per_grid_point%data1D(i) == 8) then
243              psigold(coord(1), 2, coord(2), 1, coord(3), 1) = wvl_band(iCoeff + 0)
244              psigold(coord(1), 1, coord(2), 2, coord(3), 1) = wvl_band(iCoeff + 1)
245              psigold(coord(1), 2, coord(2), 2, coord(3), 1) = wvl_band(iCoeff + 2)
246              psigold(coord(1), 1, coord(2), 1, coord(3), 2) = wvl_band(iCoeff + 3)
247              psigold(coord(1), 2, coord(2), 1, coord(3), 2) = wvl_band(iCoeff + 4)
248              psigold(coord(1), 1, coord(2), 2, coord(3), 2) = wvl_band(iCoeff + 5)
249              psigold(coord(1), 2, coord(2), 2, coord(3), 2) = wvl_band(iCoeff + 6)
250              iCoeff = iCoeff + 7
251            end if
252          end do
253 
254          call reformatonewave(displ, wfs%ks%lzd%Glr%wfd, wvl%atoms, hgrid_old(1), &
255 &         hgrid_old(2), hgrid_old(3), n_old(1), n_old(2), &
256 &         n_old(3), xcart_old, psigold, &
257 &         wvl%h(1), wvl%h(2), wvl%h(3), &
258 &         wvl%Glr%d%n1, wvl%Glr%d%n2, &
259 &         wvl%Glr%d%n3, xcart, psifscf, &
260 &         wfs%ks%psi(bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + 1: &
261 &         bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + bandSize))
262        else
263 !        No reformating, just reorder
264          iCoeff = 1
265          iFine = bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + wfs%ks%lzd%Glr%wfd%nvctr_c + 1
266          do iCoarse = 1, hdr0%nwvlarr(1), 1
267            wfs%ks%psi(bandSize * (iBand - me * wfs%ks%orbs%norbp - 1) + iCoarse) = &
268 &           wvl_band(iCoeff)
269            iCoeff  = iCoeff  + 1
270            if (basis_folder%number_of_coefficients_per_grid_point%data1D(iCoarse) == 8) then
271              wfs%ks%psi(iFine:iFine + 6) = wvl_band(iCoeff:iCoeff + 6)
272              iFine  = iFine  + 7
273              iCoeff = iCoeff + 7
274            end if
275          end do
276        end if
277      end do
278 
279 !    We read wfs%ks%eval (TODO maybe removed later).
280      electrons_folder%eigenvalues%data1D => wfs%ks%orbs%eval
281      call etsf_io_electrons_get(wff%unwff, electrons_folder, lstat, error)
282      ETSF_CHECK_ERROR(lstat, error)
283 
284 !    Deallocate temporary arrays
285      ABI_DEALLOCATE(wvl_band)
286      ABI_DEALLOCATE(basis_folder%number_of_coefficients_per_grid_point%data1D)
287      if (reformat) then
288        ABI_DEALLOCATE(basis_folder%coordinates_of_basis_grid_points%data2D)
289        ABI_DEALLOCATE(psigold)
290        ABI_DEALLOCATE(psifscf)
291      end if
292 #endif
293 
294    else
295      write(message,'(4a,i0,a)') ch10,&
296 &     '  wff%iomode argument is wrong,', ch10, &
297 &     '  awaited values are -1, 0 (or 3 if netcdf/etsf_io is available) but value = ', wff%iomode, '.'
298      MSG_BUG(message)
299    end if
300  else
301    call readmywaves(me, "wavefunctions", WF_FORMAT_NONE, wfs%ks%orbs, &
302 &   wvl%Glr%d%n1, wvl%Glr%d%n2, wvl%Glr%d%n3, &
303 &   wvl%h(1), wvl%h(2), wvl%h(3), wvl%atoms, &
304 &   xcart_old, xcart, wfs%ks%lzd%Glr%wfd, wfs%ks%psi)
305  end if
306 
307  ABI_DEALLOCATE(xcart)
308  ABI_DEALLOCATE(xcart_old)
309 #else
310  BIGDFT_NOTENABLED_ERROR()
311  if (.false.) write(std_out,*) option,dtset%nstep,hdr0%ecut,hdr%ecut,mpi_enreg%nproc,wff%me,&
312 & wfs%ks,wvl%h(1),rprimd(1,1),xred(1,1)
313 #endif
314 
315 end subroutine wvl_read

ABINIT/wvl_write [ Functions ]

[ Top ] [ 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

350 subroutine wvl_write(dtset, eigen, mpi_enreg, option, rprimd, wff, wfs, wvl, xred)
351 
352 #if defined HAVE_BIGDFT
353   use BigDFT_API, only : writeonewave,writemywaves,WF_FORMAT_NONE
354 #endif
355 
356 !This section has been created automatically by the script Abilint (TD).
357 !Do not modify the following lines by hand.
358 #undef ABI_FUNC
359 #define ABI_FUNC 'wvl_write'
360 !End of the abilint section
361 
362   implicit none
363 
364 !Arguments -------------------------------
365   !scalars
366   integer, intent(in)                       :: option
367   type(dataset_type), intent(in)            :: dtset
368   type(MPI_type), intent(in)                :: mpi_enreg
369   type(wffile_type),intent(in)              :: wff
370   type(wvl_wf_type), intent(in)             :: wfs
371   type(wvl_internal_type), intent(in)       :: wvl
372   !arrays
373   real(dp), intent(in), target              :: eigen(dtset%mband)
374   real(dp), intent(in)                      :: rprimd(3, 3)
375   real(dp), intent(in)                      :: xred(3, dtset%natom)
376 
377 !Local variables-------------------------------
378 #if defined HAVE_BIGDFT
379   character(len = 500)  :: message
380   integer               :: comm,me
381   integer               :: iorb
382   integer               :: iseg, nseg, ipsi, npsi
383   real(dp), allocatable :: xcart(:,:)
384 #if defined HAVE_ETSF_IO
385   integer               :: i, i0, i1, i2, i3, jj, j0, j1, ii
386   integer               :: iGrid, iCoeff, iCoarse, iFine
387   integer               :: iband
388   integer, allocatable  :: coeff_map(:,:,:)
389   integer, allocatable, target :: wvl_coord(:,:), wvl_ncoeff(:)
390   real(dp), allocatable, target :: wvl_band(:)
391   type(etsf_basisdata) :: basis_folder
392   type(etsf_main)      :: main_folder
393   type(etsf_electrons) :: electrons_folder
394   logical              :: lstat
395   type(etsf_io_low_error) :: error
396 #endif
397 #endif
398 
399 ! *********************************************************************
400 
401 #if defined HAVE_BIGDFT
402 
403  if (abs(option) /= 2) then
404    write(message,'(a,a,a,i0,a)')&
405 &   '  Option argument is wrong,', ch10, &
406 &   '  awaited values are -2 or  2 but option = ', option, '.'
407    MSG_BUG(message)
408  end if
409 
410  comm=mpi_enreg%comm_wvl
411  me=xmpi_comm_rank(comm)
412 !Store xcart for each atom
413  ABI_ALLOCATE(xcart,(3, dtset%natom))
414  call xred2xcart(dtset%natom, rprimd, xcart, xred)
415 
416  write(message, '(a,a,a,a)' ) ch10,&
417 & ' wvl_write:  Write wavefunctions to file.'
418  call wrtout(std_out,message,'COLL')
419 
420  if (option > 0) then
421 !  Write in the ABINIT way.
422    if (wff%iomode == IO_MODE_FORTRAN .or. (wff%iomode == IO_MODE_FORTRAN_MASTER .and. wff%master==wff%me)) then
423      iseg = wfs%ks%lzd%Glr%wfd%nseg_c
424      nseg = wfs%ks%lzd%Glr%wfd%nseg_c + wfs%ks%lzd%Glr%wfd%nseg_f
425      ipsi = wfs%ks%lzd%Glr%wfd%nvctr_c
426      npsi = wfs%ks%lzd%Glr%wfd%nvctr_c + 7 * wfs%ks%lzd%Glr%wfd%nvctr_f
427      do iorb = 1, dtset%mband
428        call writeonewave(wff%unwff, .false., iorb, wvl%Glr%d%n1, &
429 &       wvl%Glr%d%n2, wvl%Glr%d%n3, &
430 &       wvl%h(1), wvl%h(2), wvl%h(3), dtset%natom, &
431 &       xcart, wfs%ks%lzd%Glr%wfd%nseg_c, wfs%ks%lzd%Glr%wfd%nvctr_c, &
432 &       wfs%ks%lzd%Glr%wfd%keygloc(:,1:iseg), wfs%ks%lzd%Glr%wfd%keyvloc(1:iseg), wfs%ks%lzd%Glr%wfd%nseg_f, &
433 &       wfs%ks%lzd%Glr%wfd%nvctr_f, wfs%ks%lzd%Glr%wfd%keygloc(:, iseg + 1:nseg), &
434 &       wfs%ks%lzd%Glr%wfd%keyvloc(iseg + 1:nseg), &
435 &       wfs%ks%psi(npsi * (iorb - me * wfs%ks%orbs%norbp - 1) + 1: &
436 &       npsi * (iorb - me * wfs%ks%orbs%norbp - 1) + ipsi), &
437 &       wfs%ks%psi(npsi * (iorb - me * wfs%ks%orbs%norbp - 1) + ipsi + 1: &
438 &       npsi * (iorb - me * wfs%ks%orbs%norbp - 1) + npsi), &
439 &       wfs%ks%orbs%eval(iorb))
440      end do
441 #if defined HAVE_ETSF_IO
442    else if (wff%iomode == IO_MODE_ETSF) then
443 !    Write a NetCDF file.
444 !    coordinates_of_grid_points is used to store the geometric
445 !    position of coefficients of wavelets i, as integer in
446 !    wvl%dpSize(:) dimensions.
447 !    coefficients_of_wavefunctions is used to store the psi values for
448 !    each wavelet.
449 
450 !    We write the basis set.
451 !    =======================
452      ABI_ALLOCATE(wvl_coord,(3, wfs%ks%lzd%Glr%wfd%nvctr_c))
453      ABI_ALLOCATE(wvl_ncoeff,(wfs%ks%lzd%Glr%wfd%nvctr_c))
454 !    Will store the grid index for a given geometric point
455      ABI_ALLOCATE(coeff_map,(0:wvl%Glr%d%n1,0:wvl%Glr%d%n2,0:wvl%Glr%d%n3))
456 !    coarse part
457      iGrid = 0
458      coeff_map = 0
459      do iseg = 1, wfs%ks%lzd%Glr%wfd%nseg_c
460        jj = wfs%ks%lzd%Glr%wfd%keyvloc(iseg)
461        j0 = wfs%ks%lzd%Glr%wfd%keygloc(1, iseg)
462        j1 = wfs%ks%lzd%Glr%wfd%keygloc(2, iseg)
463        ii = j0 - 1
464        i3 = ii / ((wvl%Glr%d%n1 + 1) * &
465 &       (wvl%Glr%d%n2 + 1))
466        ii = ii - i3 * (wvl%Glr%d%n1 + 1) * &
467 &       (wvl%Glr%d%n2 + 1)
468        i2 = ii / (wvl%Glr%d%n1 + 1)
469        i0 = ii - i2 * (wvl%Glr%d%n1 + 1)
470        i1 = i0 + j1 - j0
471        do i = i0, i1
472          iGrid = iGrid + 1
473          coeff_map(i, i2, i3) = iGrid
474          wvl_coord(:, iGrid) = (/ i * 2, i2 * 2, i3 * 2 /)
475          wvl_ncoeff(iGrid) = 1
476        end do
477      end do
478 !    fine part
479      do iseg = 1, wfs%ks%lzd%Glr%wfd%nseg_f
480        jj = wfs%ks%lzd%Glr%wfd%keyvloc(wfs%ks%lzd%Glr%wfd%nseg_c + iseg)
481        j0 = wfs%ks%lzd%Glr%wfd%keygloc(1, wfs%ks%lzd%Glr%wfd%nseg_c + iseg)
482        j1 = wfs%ks%lzd%Glr%wfd%keygloc(2, wfs%ks%lzd%Glr%wfd%nseg_c + iseg)
483        ii = j0 - 1
484        i3 = ii / ((wvl%Glr%d%n1 + 1) * &
485 &       (wvl%Glr%d%n2 + 1))
486        ii = ii - i3 * (wvl%Glr%d%n1 + 1) * &
487 &       (wvl%Glr%d%n2 + 1)
488        i2 = ii / (wvl%Glr%d%n1 + 1)
489        i0 = ii - i2 * (wvl%Glr%d%n1 + 1)
490        i1 = i0 + j1 - j0
491        do i = i0, i1
492          iGrid = coeff_map(i , i2 , i3)
493          wvl_ncoeff(iGrid) = wvl_ncoeff(iGrid) + 7
494        end do
495      end do
496      basis_folder%coordinates_of_basis_grid_points%data2D => wvl_coord
497      basis_folder%number_of_coefficients_per_grid_point%data1D => wvl_ncoeff
498 
499      call etsf_io_basisdata_put(wff%unwff, basis_folder, lstat, error)
500 
501      ETSF_CHECK_ERROR(lstat, error)
502 
503 !    We write the wavefuctions.
504 !    ==========================
505 !    We write the wavefunctions band per band.
506      ABI_ALLOCATE(wvl_band,(wfs%ks%lzd%Glr%wfd%nvctr_c + 7 * wfs%ks%lzd%Glr%wfd%nvctr_f))
507 !    We reorder the wfs%ks%psi(:,iband) into wvl_band
508      do iband = 1, dtset%mband, 1
509 !      iCoeff is the index of the coefficient we are writing in wvl_band
510        iCoeff  = 1
511 !      iCoarse is the index of the coarse part in wfs%ks%psi (=iGrid)
512        iCoarse = 1
513 !      iFine is the index of the fine part in wfs%ks%psi
514        iFine   = wfs%ks%lzd%Glr%wfd%nvctr_c + 1
515        npsi    = wfs%ks%lzd%Glr%wfd%nvctr_c + 7 * wfs%ks%lzd%Glr%wfd%nvctr_f
516        do iGrid = 1, wfs%ks%lzd%Glr%wfd%nvctr_c, 1
517          wvl_band(iCoeff) = &
518 &         wfs%ks%psi((iband - me * wfs%ks%orbs%norbp - 1) * npsi + iCoarse)
519          iCoeff  = iCoeff + 1
520          iCoarse = iCoarse + 1
521          if (wvl_ncoeff(iGrid) == 8) then
522            wvl_band(iCoeff:iCoeff + 6) = &
523 &           wfs%ks%psi((iband - me * wfs%ks%orbs%norbp - 1) * npsi + iFine: &
524 &           (iband - me * wfs%ks%orbs%norbp - 1) * npsi + iFine + 6)
525            iCoeff = iCoeff + 7
526            iFine  = iFine + 7
527          end if
528        end do
529        main_folder%wfs_coeff__state_access = iband
530        main_folder%coefficients_of_wavefunctions%data1D => wvl_band
531        call etsf_io_main_put(wff%unwff, main_folder, lstat, error)
532        ETSF_CHECK_ERROR(lstat, error)
533      end do
534 
535 !    We write the electronic informations.
536 !    =====================================
537      electrons_folder%eigenvalues%data1D => eigen
538      electrons_folder%eigenvalues__number_of_states = dtset%mband
539      call etsf_io_electrons_put(wff%unwff, electrons_folder, lstat, error)
540      ETSF_CHECK_ERROR(lstat,error)
541 
542 !    We deallocate all arrays
543      ABI_DEALLOCATE(wvl_band)
544      ABI_DEALLOCATE(wvl_coord)
545      ABI_DEALLOCATE(wvl_ncoeff)
546      ABI_DEALLOCATE(coeff_map)
547 #endif
548    else
549      write(message,'(3a,i0,a)')&
550 &     '  wff%iomode argument is wrong,', ch10, &
551 &     '  awaited values are -1, 0 (or 3 if netcdf/etsf_io is available) but value = ', wff%iomode, '.'
552      MSG_BUG(message)
553    end if
554  else
555 !  Write in the BigDFT way.
556    call  writemywaves(me, "wavefunctions", WF_FORMAT_NONE, wfs%ks%orbs, &
557 &   wvl%Glr%d%n1, wvl%Glr%d%n2, wvl%Glr%d%n3, &
558 &   wvl%h(1), wvl%h(2), wvl%h(3),wvl%atoms, &
559 &   xcart, wfs%ks%lzd%Glr%wfd, wfs%ks%psi)
560  end if
561 
562  ABI_DEALLOCATE(xcart)
563 
564 #else
565  BIGDFT_NOTENABLED_ERROR()
566  if (.false.) write(std_out,*) option,dtset%nstep,mpi_enreg%nproc,wff%me,&
567 & wfs%ks,wvl%h(1),eigen(1),rprimd(1,1),xred(1,1)
568 #endif
569 
570 end subroutine wvl_write