TABLE OF CONTENTS


ABINIT/m_fourier_interpol [ Modules ]

[ Top ] [ Modules ]

NAME

  m_fourier_interpol

FUNCTION

  This module contains routines used to perform a Fourier interpolation.
  Mainly used in PAW to interpol data from/to the coarse FFT grid from/to the fine FFT grid.

COPYRIGHT

 Copyright (C) 2018-2018 ABINIT group (FJ, MT, MG)
 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

18 #if defined HAVE_CONFIG_H
19 #include "config.h"
20 #endif
21 
22 #include "abi_common.h"
23 
24 MODULE m_fourier_interpol
25 
26  use defs_basis
27  use defs_abitypes
28  use m_abicore
29  use m_errors
30 
31  use m_fft,    only : zerosym, indirect_parallel_Fourier, fourdp
32  use m_pawfgr, only : pawfgr_type,pawfgr_destroy,indgrid
33 
34  implicit none
35 
36  private
37 
38 !public procedures.
39  public :: transgrid        ! Convert a density/potential from the coarse to the fine grid and vice versa
40  public :: fourier_interpol ! Perform a Fourier interpolation. Just a wrapper for transgrid.
41 
42 CONTAINS  !========================================================================================

m_fourier_interpol/fourier_interpol [ Functions ]

[ Top ] [ m_fourier_interpol ] [ Functions ]

NAME

 fourier_interpol

FUNCTION

  Perform a Fourier interpolation. Just a wrapper for transgrid, the table giving the correspondence 
  between the coarse and the mesh FFT grids are constructed inside the routine. This allows to specify an 
  arbitrary FFT mesh to be used for the interpolation. Besides the routine works also in 
  case of NC calculations since it does not require Pawfgr.

INPUTS

 cplex=1 if rhor[f] is real, 2 if rhor[f] is complex
 MPI_enreg<MPI_type>=Information about MPI parallelization
 nspden=number of spin-density components
 nfft_in =number of points in the input FFT box (WARNING no FFT parallelism)
 nfft_out=number of points in the output FFT box
 ngfft_in(18)=all needed information about 3D FFT, for the input grid
 ngfft_out(18) =all needed information about 3D FFT, for the output grid
 optin= 0: input density/potential is taken from rhor_in(:,nspden)
        1: input density/potential is taken from rhog_in(:)     (ispden=1)
                                             and rhor_in(:,2:4) (ispden=2,3,4)
 optout= 0: output density/potential is given in r space in rhor_out(:,nspden)
         1: output density/potential is given in r space in rhor_out(:,nspden)
                                          and in g space in rhog_out(:)
 ngfft_asked(18)=All info on the required FFT mesh.
 paral_kgb=Flag related to the kpoint-band-fft parallelism (TODO not implemented)

OUTPUT

  rhor_out(cplex*nfft_out,nspden)=output density/potential in r space on the required FFT mesh.
  if optout=1:
   rhog_out(2,nfftc)=Fourier transform of output density/potential on the coarse grid

PARENTS

      m_dvdb,m_ioarr,m_qparticles

CHILDREN

      indgrid,pawfgr_destroy,transgrid

SOURCE

438 subroutine fourier_interpol(cplex,nspden,optin,optout,nfft_in,ngfft_in,nfft_out,ngfft_out,&
439 & paral_kgb,MPI_enreg,rhor_in,rhor_out,rhog_in,rhog_out)
440 
441 
442 !This section has been created automatically by the script Abilint (TD).
443 !Do not modify the following lines by hand.
444 #undef ABI_FUNC
445 #define ABI_FUNC 'fourier_interpol'
446 !End of the abilint section
447 
448  implicit none
449 
450 !Arguments ------------------------------------
451 !scalars
452  integer,intent(in) :: cplex,nspden,optin,optout
453  integer,intent(in) :: nfft_in,nfft_out,paral_kgb
454  type(MPI_type),intent(in) :: MPI_enreg
455 !arrays
456  integer,intent(in) :: ngfft_in(18),ngfft_out(18)
457  real(dp),intent(inout) :: rhor_in(cplex*nfft_in,nspden)
458  real(dp),intent(inout) :: rhog_in(2,nfft_in)
459  real(dp),intent(out) :: rhor_out(cplex*nfft_out,nspden)
460  real(dp),intent(out) :: rhog_out(2,nfft_out)
461 
462 !Local variables ---------------------------------------
463 !scalars
464  integer :: nfftf,nfftc,nfftc_tot,nfftf_tot,optgrid
465  logical :: ltest
466  type(Pawfgr_type) :: Pawfgr
467 !arrays
468  integer :: ngfftc(18),ngfftf(18)
469 
470 ! *************************************************************************
471 
472 !=== FFT parallelism not implemented ===
473  ABI_CHECK(paral_kgb==0,'paral_kgb/=0 not implemented')
474 
475  ltest= ALL( ngfft_in(7:) == ngfft_out(7:) )
476  ABI_CHECK(ltest,'ngfftf_in(7:18)/=ngfftf_out(7:18)')
477 
478 !================================
479 !=== Which one is the coarse? ===
480 !================================
481  if (nfft_out>=nfft_in) then 
482    ! From coarse to fine grid. If meshes are equivalent, call transgrid anyway because of optout, optin.
483    nfftf    =nfft_out 
484    ngfftf(:)=ngfft_out(:)
485    nfftf_tot =PRODUCT(ngfft_out(1:3))
486 
487    nfftc    =nfft_in 
488    ngfftc(:)=ngfft_in(:)
489    nfftc_tot =PRODUCT(ngfft_in (1:3))
490 
491    Pawfgr%usefinegrid=1 
492    if (ALL(ngfft_in(1:3)==ngfft_out(1:3))) Pawfgr%usefinegrid=0 
493    optgrid=1
494 
495  else 
496    ! From fine towards coarse.
497    nfftf    =nfft_in 
498    ngfftf(:)=ngfft_in(:)
499    nfftf_tot =PRODUCT(ngfft_in(1:3))
500 
501    nfftc    =nfft_out 
502    ngfftc(:)=ngfft_out(:)
503    nfftc_tot =PRODUCT(ngfft_out (1:3))
504    Pawfgr%usefinegrid=1 
505    optgrid=-1
506  end if
507 
508  ABI_MALLOC(Pawfgr%coatofin,(nfftc_tot))
509  ABI_MALLOC(Pawfgr%fintocoa,(nfftf_tot))
510 
511  call indgrid(Pawfgr%coatofin,Pawfgr%fintocoa,nfftc_tot,nfftf_tot,ngfftc,ngfftf)
512 
513  Pawfgr%mgfft =MAXVAL (ngfftf(1:3))
514  Pawfgr%nfft  =PRODUCT(ngfftf(1:3)) ! no FFT parallelism!
515  Pawfgr%ngfft =ngfftf
516 
517  Pawfgr%mgfftc =MAXVAL (ngfftc(1:3))
518  Pawfgr%nfftc  =PRODUCT(ngfftc(1:3)) ! no FFT parallelism!
519  Pawfgr%ngfftc =ngfftc
520 
521  if (optgrid==1) then 
522    call transgrid(cplex,MPI_enreg,nspden,optgrid,optin,optout,paral_kgb,Pawfgr,rhog_in ,rhog_out,rhor_in ,rhor_out)
523  else
524    call transgrid(cplex,MPI_enreg,nspden,optgrid,optin,optout,paral_kgb,Pawfgr,rhog_out,rhog_in ,rhor_out,rhor_in)
525  end if
526 
527  call pawfgr_destroy(Pawfgr)
528  
529 end subroutine fourier_interpol

m_fourier_interpol/transgrid [ Functions ]

[ Top ] [ m_fourier_interpol ] [ Functions ]

NAME

 transgrid

FUNCTION

 Convert a given density (or potential) from the coarse to the fine rectangular grid and vice versa
 Used in PAW calculations

INPUTS

  cplex=1 if rhor[f] is real, 2 if rhor[f] is complex
  mpi_enreg=information about MPI parallelization
  nspden=number of spin-density components
  optgrid=+1 to go from the coarse grid towards the fine grid
          -1 to go from the fine grid towards the coarse grid
  optin= 0: input density/potential is taken from rhor(:,nspden)
         1: input density/potential is taken from rhog(:)     (ispden=1)
                                              and rhor(:,2:4) (ispden=2,3,4)
  optout= 0: output density/potential is given in r space in rhor(:,nspden)
          1: output density/potential is given in r space in rhor(:,nspden)
                                           and in g space in rhog(:)
  pawfgr <type(paw_fgr_type)>=fine rectangular grid parameters
    %nfftc=number of points in the coarse FFT box
    %nfft =number of points in the fine FFT box
    %ngfftc(18)=all needed information about 3D FFT, for the coarse grid
    %ngfft(18) =all needed information about 3D FFT, for the fine grid
    %coatofin(nfftc)=Index of the points of the coarse grid on the fine grid
    %fintocoa(nfft) =Index of the points of the fine grid on the coarse grid
    %usefinegrid= 1 if a fine FFT grid is used (0 otherwise)
  if optgrid=+1 and optin=1:
    rhog(2,nfftc)=Fourier transform of input density/potential on the coarse grid
  if optgrid=-1 and optin=1:
    rhogf(2,nfftf)=Fourier transform of input density/potential on the fine grid
  if optgrid=+1
    rhor(cplex*nfftc,nspden)=input density/potential in r space on the coarse grid
  if optgrid=-1:
    rhorf(cplex*nfftf,nspden)=input density/potential in r space on the fine grid

OUTPUT

  if optgrid=-1 and optout=1:
    rhog(2,nfftc)=Fourier transform of output density/potential on the coarse grid
  if optgrid=+1 and optout=1:
    rhogf(2,nfftf)=Fourier transform of output density/potential on the fine grid
  if optgrid=-1
    rhor(cplex*nfftc,nspden)=output density/potential in r space on the coarse grid
  if optgrid=+1:
    rhorf(cplex*nfftf,nspden)=output density/potential in r space on the fine grid

PARENTS

      dfpt_looppert,energy,fourier_interpol,getgh1c,gstate,ks_ddiago,m_io_kss
      pawmkrho,respfn,scfcv,vtorho,vtorhorec

CHILDREN

      fourdp,indirect_parallel_fourier,zerosym

SOURCE

103 subroutine transgrid(cplex,mpi_enreg,nspden,optgrid,optin,optout,paral_kgb,pawfgr,rhog,rhogf,rhor,rhorf)
104 
105 
106 !This section has been created automatically by the script Abilint (TD).
107 !Do not modify the following lines by hand.
108 #undef ABI_FUNC
109 #define ABI_FUNC 'transgrid'
110 !End of the abilint section
111 
112  implicit none
113 
114 !Arguments ---------------------------------------------
115 !scalars
116  integer,intent(in) :: cplex,nspden,optgrid,optin,optout,paral_kgb
117  type(MPI_type),intent(in) :: mpi_enreg
118  type(pawfgr_type),intent(in) :: pawfgr
119 !arrays
120  real(dp),intent(inout) :: rhog(2,pawfgr%nfftc),rhogf(2,pawfgr%nfft)
121  real(dp),intent(inout) :: rhor(cplex*pawfgr%nfftc,nspden),rhorf(cplex*pawfgr%nfft,nspden)
122 
123 !Local variables ---------------------------------------
124 !scalars
125  integer :: i1,ispden,nfftc,nfftctot,nfftf,nfftftot
126  character(len=500) :: msg
127 !arrays
128  integer :: ngfftc(18),ngfftf(18)
129  real(dp),allocatable :: vectg(:,:),work(:,:),workfft(:)
130 
131 ! *************************************************************************
132 
133  DBG_ENTER("COLL")
134 
135 !Tests
136  if(pawfgr%nfft<pawfgr%nfftc) then
137    write(msg,'(a,2(i0,1x))')' nfft (fine grid) must be >= nfft (coarse grid) while: ',pawfgr%nfft, pawfgr%nfftc
138    MSG_ERROR(msg)
139  end if
140 
141 !Store FFT dimensions
142  nfftc=pawfgr%nfftc;ngfftc(:)=pawfgr%ngfftc(:);nfftctot=ngfftc(1)*ngfftc(2)*ngfftc(3)
143  nfftf=pawfgr%nfft ;ngfftf(:)=pawfgr%ngfft (:);nfftftot=ngfftf(1)*ngfftf(2)*ngfftf(3)
144 
145 !If no fine FFT grid is used, this is only a simple transfer
146  if (pawfgr%usefinegrid==0) then
147    if (optgrid==1) then
148      rhorf=rhor
149      if (optout==1.and.optin==1) rhogf=rhog
150      if (optout==1.and.optin/=1) then
151        ABI_ALLOCATE(workfft,(cplex*nfftc))
152        workfft(:)=rhor(:,1)
153        call fourdp(cplex,rhogf,workfft,-1,mpi_enreg,nfftc,ngfftc,paral_kgb,0)
154        ABI_DEALLOCATE(workfft)
155      end if
156    end if
157    if (optgrid==-1) then
158      rhor=rhorf
159      if (optout==1.and.optin==1) rhog=rhogf
160      if (optout==1.and.optin/=1) then
161        ABI_ALLOCATE(workfft,(cplex*nfftc))
162        workfft(:)=rhorf(:,1)
163        call fourdp(cplex,rhog,workfft,-1,mpi_enreg,nfftc,ngfftc,paral_kgb,0)
164        ABI_DEALLOCATE(workfft)
165      end if
166    end if
167    return
168  end if
169 
170 !====== FROM THE COARSE GRID TOWARDS THE FINE GRID =============
171 !===============================================================
172 !Calculate the FT of rhor to have it in the g space on the coarse grid
173 !Transfer the FT of rhor on the coarse grid towards the fine grid
174 !Then calculate the FT back to get rhorf on the fine grid
175  if (optgrid==1) then
176 
177    ABI_ALLOCATE(work,(2,nfftc))
178 
179 !  First spin component
180 !  --------------------------------------------------------------
181    if (optout==0) then
182 !    if optout=0, rhog on the fine grid is temporary (in vectg)
183      ABI_ALLOCATE(vectg,(2,nfftf))
184      vectg(:,:)=zero
185      if (optin==1) then
186        call zerosym(rhog,2,ngfftc(1),ngfftc(2),ngfftc(3),&
187 &       comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
188        if(mpi_enreg%nproc_fft > 1 .and. mpi_enreg%paral_kgb==1) then
189          call indirect_parallel_Fourier&
190 &         (pawfgr%coatofin,vectg,mpi_enreg,ngfftf,ngfftc,nfftf,nfftc,paral_kgb,rhog,nfftctot)
191        else
192          do i1=1,nfftc
193            vectg(:,pawfgr%coatofin(i1))=rhog(:,i1)
194          end do
195        end if
196      else
197        ABI_ALLOCATE(workfft,(cplex*nfftc))
198        workfft(:)=rhor(:,1)
199        call fourdp(cplex,work,workfft,-1,mpi_enreg,nfftc,ngfftc,paral_kgb,0)
200        ABI_DEALLOCATE(workfft)
201        call zerosym(work,2,ngfftc(1),ngfftc(2),ngfftc(3),&
202 &       comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
203        if(mpi_enreg%nproc_fft > 1 .and. mpi_enreg%paral_kgb==1) then
204          call indirect_parallel_Fourier&
205 &         (pawfgr%coatofin,vectg,mpi_enreg,ngfftf,ngfftc,nfftf,nfftc,paral_kgb,work,nfftctot)
206        else
207          do i1=1,nfftc
208            vectg(:,pawfgr%coatofin(i1))=work(:,i1)
209          end do
210        end if
211      end if
212 !    call zerosym(vectg,2,ngfftf(1),ngfftf(2),ngfftf(3),&
213 !    &        comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
214      ABI_ALLOCATE(workfft,(cplex*nfftf))
215      call fourdp(cplex,vectg,workfft,1,mpi_enreg,nfftf,ngfftf,paral_kgb,0)
216      rhorf(:,1)=workfft(:)
217      ABI_DEALLOCATE(workfft)
218      ABI_DEALLOCATE(vectg)
219    else
220 !    if optout=1, rhog on the fine grid is saved
221      call zerosym(rhog,2,ngfftc(1),ngfftc(2),ngfftc(3),&
222 &     comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
223      rhogf(:,:)=zero
224      if (optin==1) then
225        if(mpi_enreg%nproc_fft > 1 .and. mpi_enreg%paral_kgb==1) then
226          call indirect_parallel_Fourier&
227 &         (pawfgr%coatofin,rhogf,mpi_enreg,ngfftf,ngfftc,nfftf,nfftc,paral_kgb,rhog,nfftctot)
228        else
229          do i1=1,nfftc
230            rhogf(:,pawfgr%coatofin(i1))=rhog(:,i1)
231          end do
232        end if
233      else
234        ABI_ALLOCATE(workfft,(cplex*nfftc))
235        workfft(:)=rhor(:,1)
236        call fourdp(cplex,work,workfft,-1,mpi_enreg,nfftc,ngfftc,paral_kgb,0)
237        ABI_DEALLOCATE(workfft)
238        call zerosym(work,2,ngfftc(1),ngfftc(2),ngfftc(3),&
239 &       comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
240        if(mpi_enreg%nproc_fft > 1 .and. mpi_enreg%paral_kgb==1) then
241          call indirect_parallel_Fourier&
242 &         (pawfgr%coatofin,rhogf,mpi_enreg,ngfftf,ngfftc,nfftf,nfftc,paral_kgb,work,nfftctot)
243        else
244          do i1=1,nfftc
245            rhogf(:,pawfgr%coatofin(i1))=work(:,i1)
246          end do
247        end if
248      end if
249 !    call zerosym(rhogf,2,ngfftf(1),ngfftf(2),ngfftf(3),&
250 !    &        comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
251      ABI_ALLOCATE(workfft,(cplex*nfftf))
252      call fourdp(cplex,rhogf,workfft,1,mpi_enreg,nfftf,ngfftf,paral_kgb,0)
253      rhorf(:,1)=workfft(:)
254      ABI_DEALLOCATE(workfft)
255    end if
256 
257 !  Additional spin components
258 !  ----------------------------------------------------
259    if (nspden>=2) then
260      ABI_ALLOCATE(vectg,(2,nfftf))
261      do ispden=2,nspden
262        vectg(:,:)=zero
263        ABI_ALLOCATE(workfft,(cplex*nfftc))
264        workfft(:)=rhor(:,ispden)
265        call fourdp(cplex,work,workfft,-1,mpi_enreg,nfftc,ngfftc,paral_kgb,0)
266        ABI_DEALLOCATE(workfft)
267        call zerosym(work,2,ngfftc(1),ngfftc(2),ngfftc(3),&
268 &       comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
269        if(mpi_enreg%nproc_fft > 1 .and. mpi_enreg%paral_kgb==1) then
270          call indirect_parallel_Fourier&
271 &         (pawfgr%coatofin,vectg,mpi_enreg,ngfftf,ngfftc,nfftf,nfftc,paral_kgb,work,nfftctot)
272        else
273          do i1=1,nfftc
274            vectg(:,pawfgr%coatofin(i1))=work(:,i1)
275          end do
276        end if
277 !      call zerosym(vectg,2,ngfftf(1),ngfftf(2),ngfftf(3),&
278 !      &          comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
279        ABI_ALLOCATE(workfft,(cplex*nfftf))
280        call fourdp(cplex,vectg,workfft,1,mpi_enreg,nfftf,ngfftf,paral_kgb,0)
281        rhorf(:,ispden)=workfft(:)
282        ABI_DEALLOCATE(workfft)
283      end do
284      ABI_DEALLOCATE(vectg)
285    end if
286 
287    ABI_DEALLOCATE(work)
288 
289 
290 !  ====== FROM THE FINE GRID TOWARDS THE COARSE GRID =============
291 !  ==============================================================
292 !  Calculate the FT of rhorf to have it in the g space on the fine grid
293 !  Transfer the FT of rhorf on the fine grid towards the coarse grid
294 !  Then calculate the FT back to get rhor on the coarse grid
295  else if (optgrid==-1) then
296 
297    ABI_ALLOCATE(work,(2,nfftf))
298 
299 !  First spin component
300 !  --------------------------------------------------------------
301    if (optout==0) then
302 !    if optout=0, rhog on the fine grid is temporary (in vectg)
303      ABI_ALLOCATE(vectg,(2,nfftc))
304      vectg(:,:)=zero
305      if (optin==1) then
306        do i1=1,nfftf
307          if (pawfgr%fintocoa(i1)/=0) vectg(:,pawfgr%fintocoa(i1))=rhogf(:,i1)
308        end do
309      else
310        ABI_ALLOCATE(workfft,(cplex*nfftf))
311        workfft(:)=rhorf(:,1)
312        call fourdp(cplex,work,workfft,-1,mpi_enreg,nfftf,ngfftf,paral_kgb,0)
313        ABI_DEALLOCATE(workfft)
314        if(mpi_enreg%nproc_fft > 1 .and. mpi_enreg%paral_kgb==1) then
315          call indirect_parallel_Fourier&
316 &         (pawfgr%fintocoa,vectg,mpi_enreg,ngfftc,ngfftf,nfftc,nfftf,paral_kgb,work,nfftftot)
317        else
318          do i1=1,nfftf
319            if (pawfgr%fintocoa(i1)/=0) vectg(:,pawfgr%fintocoa(i1))=work(:,i1)
320          end do
321        end if
322      end if
323      call zerosym(vectg,2,ngfftc(1),ngfftc(2),ngfftc(3),&
324 &     comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
325      ABI_ALLOCATE(workfft,(cplex*nfftc))
326      call fourdp(cplex,vectg,workfft,1,mpi_enreg,nfftc,ngfftc,paral_kgb,0)
327      rhor(:,1)=workfft(:)
328      ABI_DEALLOCATE(workfft)
329      ABI_DEALLOCATE(vectg)
330    else
331 !    if optout=1, rhog on the fine grid is saved
332      rhog(:,:)=zero
333      if (optin==1) then
334        do i1=1,nfftf
335          if (pawfgr%fintocoa(i1)/=0) rhog(:,pawfgr%fintocoa(i1))=rhogf(:,i1)
336        end do
337      else
338        ABI_ALLOCATE(workfft,(cplex*nfftf))
339        workfft(:)=rhorf(:,1)
340        call fourdp(cplex,work,workfft,-1,mpi_enreg,nfftf,ngfftf,paral_kgb,0)
341        ABI_DEALLOCATE(workfft)
342        if(mpi_enreg%nproc_fft > 1 .and. mpi_enreg%paral_kgb==1) then
343          call indirect_parallel_Fourier&
344 &         (pawfgr%fintocoa,rhog,mpi_enreg,ngfftc,ngfftf,nfftc,nfftf,paral_kgb,work,nfftftot)
345        else
346          do i1=1,nfftf
347            if (pawfgr%fintocoa(i1)/=0) rhog(:,pawfgr%fintocoa(i1))=work(:,i1)
348          end do
349        end if
350      end if
351      call zerosym(rhog,2,ngfftc(1),ngfftc(2),ngfftc(3),&
352 &     comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
353      ABI_ALLOCATE(workfft,(cplex*nfftc))
354      call fourdp(cplex,rhog,workfft,1,mpi_enreg,nfftc,ngfftc,paral_kgb,0)
355      rhor(:,1)=workfft(:)
356      ABI_DEALLOCATE(workfft)
357    end if
358 
359 !  Additional spin components
360 !  ----------------------------------------------------
361    if (nspden>=2) then
362      ABI_ALLOCATE(vectg,(2,nfftc))
363      do ispden=2,nspden
364        vectg(:,:)=zero
365        ABI_ALLOCATE(workfft,(cplex*nfftf))
366        workfft(:)=rhorf(:,ispden)
367        call fourdp(cplex,work,workfft,-1,mpi_enreg,nfftf,ngfftf,paral_kgb,0)
368        ABI_DEALLOCATE(workfft)
369        if(mpi_enreg%nproc_fft > 1 .and. mpi_enreg%paral_kgb==1) then
370          call indirect_parallel_Fourier&
371 &         (pawfgr%fintocoa,vectg,mpi_enreg,ngfftc,ngfftf,nfftc,nfftf,paral_kgb,work,nfftftot)
372        else
373          do i1=1,nfftf
374            if (pawfgr%fintocoa(i1)/=0) vectg(:,pawfgr%fintocoa(i1))=work(:,i1)
375          end do
376        end if
377        call zerosym(vectg,2,ngfftc(1),ngfftc(2),ngfftc(3),&
378 &       comm_fft=mpi_enreg%comm_fft,distribfft=mpi_enreg%distribfft)
379        ABI_ALLOCATE(workfft,(cplex*nfftc))
380        call fourdp(cplex,vectg,workfft,1,mpi_enreg,nfftc,ngfftc,paral_kgb,0)
381        rhor(:,ispden)=workfft(:)
382        ABI_DEALLOCATE(workfft)
383      end do
384      ABI_DEALLOCATE(vectg)
385    end if
386 
387    ABI_DEALLOCATE(work)
388 
389  end if
390 
391  DBG_EXIT("COLL")
392 
393 end subroutine transgrid