TABLE OF CONTENTS


ABINIT/m_distribfft [ Modules ]

[ Top ] [ Modules ]

NAME

  m_distribfft

FUNCTION

  This module provides the definition of the different arrays
  used for FFT parallelization with MPI and n2 plane sharing

COPYRIGHT

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

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 #include "abi_common.h"
27 
28 MODULE m_distribfft
29 
30  use defs_basis
31  use m_abicore
32  use m_errors
33 
34  implicit none
35 
36  private

m_distribfft/copy_distribfft [ Functions ]

[ Top ] [ m_distribfft ] [ Functions ]

NAME

  copy_distribfft

FUNCTION

  Cleans-up the mpi information for FFT distribution
  (mostly deallocate parts distribfft(:) ).

SIDE EFFECTS

  mpi_enreg=information about MPI parallelization

PARENTS

      m_mpinfo

CHILDREN

SOURCE

495 subroutine copy_distribfft(distribfft_src, distribfft_dst)
496 
497 
498 !This section has been created automatically by the script Abilint (TD).
499 !Do not modify the following lines by hand.
500 #undef ABI_FUNC
501 #define ABI_FUNC 'copy_distribfft'
502 !End of the abilint section
503 
504  implicit none
505 
506 !Arguments ------------------------------------
507  type(distribfft_type),intent(in)   :: distribfft_src
508  type(distribfft_type),intent(out) :: distribfft_dst
509 
510 ! ***********************************************************************
511 
512  DBG_ENTER("COLL")
513 
514  distribfft_dst%n2_coarse=  distribfft_src%n2_coarse
515  distribfft_dst%n2_fine  =  distribfft_src%n2_fine
516 
517  if (allocated(distribfft_src%tab_fftwf2_distrib)) then
518    ABI_ALLOCATE(distribfft_dst%tab_fftwf2_distrib,(size(distribfft_src%tab_fftwf2_distrib)))
519    distribfft_dst%tab_fftwf2_distrib=distribfft_src%tab_fftwf2_distrib
520  end if
521 
522  if (allocated(distribfft_src%tab_fftdp2_distrib)) then
523   ABI_ALLOCATE(distribfft_dst%tab_fftdp2_distrib,(size(distribfft_src%tab_fftdp2_distrib)))
524   distribfft_dst%tab_fftdp2_distrib=distribfft_src%tab_fftdp2_distrib
525  end if
526 
527  if (allocated(distribfft_src%tab_fftdp3_distrib)) then
528   ABI_ALLOCATE(distribfft_dst%tab_fftdp3_distrib,(size(distribfft_src%tab_fftdp3_distrib)))
529   distribfft_dst%tab_fftdp3_distrib=distribfft_src%tab_fftdp3_distrib
530  end if
531 
532  if (allocated(distribfft_src%tab_fftwf2dg_distrib)) then
533    ABI_ALLOCATE(distribfft_dst%tab_fftwf2dg_distrib,(size(distribfft_src%tab_fftwf2dg_distrib)))
534    distribfft_dst%tab_fftwf2dg_distrib=distribfft_src%tab_fftwf2dg_distrib
535  end if
536 
537  if (allocated(distribfft_src%tab_fftdp2dg_distrib)) then
538    ABI_ALLOCATE(distribfft_dst%tab_fftdp2dg_distrib,(size(distribfft_src%tab_fftdp2dg_distrib)))
539    distribfft_dst%tab_fftdp2dg_distrib=distribfft_src%tab_fftdp2dg_distrib
540  end if
541 
542  if (allocated(distribfft_src%tab_fftdp3dg_distrib)) then
543    ABI_ALLOCATE(distribfft_dst%tab_fftdp3dg_distrib,(size(distribfft_src%tab_fftdp3dg_distrib)))
544    distribfft_dst%tab_fftdp3dg_distrib=distribfft_src%tab_fftdp3dg_distrib
545  end if
546 
547  if (allocated(distribfft_src%tab_fftwf2_local)) then
548    ABI_ALLOCATE(distribfft_dst%tab_fftwf2_local,(size(distribfft_src%tab_fftwf2_local)))
549    distribfft_dst%tab_fftwf2_local=distribfft_src%tab_fftwf2_local
550  end if
551 
552  if (allocated(distribfft_src%tab_fftdp2_local)) then
553    ABI_ALLOCATE(distribfft_dst%tab_fftdp2_local,(size(distribfft_src%tab_fftdp2_local)))
554    distribfft_dst%tab_fftdp2_local=distribfft_src%tab_fftdp2_local
555  end if
556 
557  if (allocated(distribfft_src%tab_fftdp3_local)) then
558    ABI_ALLOCATE(distribfft_dst%tab_fftdp3_local,(size(distribfft_src%tab_fftdp3_local)))
559    distribfft_dst%tab_fftdp3_local=distribfft_src%tab_fftdp3_local
560  end if
561 
562  if (allocated(distribfft_src%tab_fftwf2dg_local)) then
563    ABI_ALLOCATE(distribfft_dst%tab_fftwf2dg_local,(size(distribfft_src%tab_fftwf2dg_local)))
564    distribfft_dst%tab_fftwf2dg_local=distribfft_src%tab_fftwf2dg_local
565  end if
566 
567  if (allocated(distribfft_src%tab_fftdp2dg_local)) then
568    ABI_ALLOCATE(distribfft_dst%tab_fftdp2dg_local,(size(distribfft_src%tab_fftdp2dg_local)))
569    distribfft_dst%tab_fftdp2dg_local=distribfft_src%tab_fftdp2dg_local
570  end if
571 
572  if (allocated(distribfft_src%tab_fftdp3dg_local)) then
573    ABI_ALLOCATE(distribfft_dst%tab_fftdp3dg_local,(size(distribfft_src%tab_fftdp3dg_local)))
574    distribfft_dst%tab_fftdp3dg_local=distribfft_src%tab_fftdp3dg_local
575  end if
576 
577  DBG_EXIT("COLL")
578 
579 end subroutine copy_distribfft

m_distribfft/destroy_distribfft [ Functions ]

[ Top ] [ m_distribfft ] [ Functions ]

NAME

  destroy_distribfft

FUNCTION

  Cleans-up the mpi information for FFT distribution
  (mostly deallocate parts distribfft(:) ).

SIDE EFFECTS

  mpi_enreg=information about MPI parallelization

PARENTS

      atm2fft,dfpt_atm2fft,m_cut3d,m_fft,m_mpinfo,m_qparticles,m_wfd
      multipoles_fftr,pawgrnl,pawmknhat,pawmknhat_psipsi,pawsushat,scfcv
      vtorhorec

CHILDREN

SOURCE

409 subroutine destroy_distribfft(distribfft_arg)
410 
411 
412 !This section has been created automatically by the script Abilint (TD).
413 !Do not modify the following lines by hand.
414 #undef ABI_FUNC
415 #define ABI_FUNC 'destroy_distribfft'
416 !End of the abilint section
417 
418  implicit none
419 
420 !Arguments ------------------------------------
421  type(distribfft_type), intent(inout) :: distribfft_arg
422 
423 ! ***********************************************************************
424 
425  DBG_ENTER("COLL")
426 
427  distribfft_arg%n2_coarse=0
428  distribfft_arg%n2_fine  =0
429 
430  if (allocated(distribfft_arg%tab_fftwf2_distrib)) then
431    ABI_DEALLOCATE(distribfft_arg%tab_fftwf2_distrib)
432  end if
433 
434  if (allocated(distribfft_arg%tab_fftdp2_distrib)) then
435    ABI_DEALLOCATE(distribfft_arg%tab_fftdp2_distrib)
436  end if
437  if (allocated(distribfft_arg%tab_fftdp3_distrib)) then
438    ABI_DEALLOCATE(distribfft_arg%tab_fftdp3_distrib)
439  end if
440  if (allocated(distribfft_arg%tab_fftwf2dg_distrib)) then
441   ABI_DEALLOCATE(distribfft_arg%tab_fftwf2dg_distrib)
442  end if
443  if (allocated(distribfft_arg%tab_fftdp2dg_distrib)) then
444   ABI_DEALLOCATE(distribfft_arg%tab_fftdp2dg_distrib)
445  end if
446  if (allocated(distribfft_arg%tab_fftdp3dg_distrib)) then
447   ABI_DEALLOCATE(distribfft_arg%tab_fftdp3dg_distrib)
448  end if
449 
450  if (allocated(distribfft_arg%tab_fftwf2_local)) then
451   ABI_DEALLOCATE(distribfft_arg%tab_fftwf2_local)
452  end if
453 
454  if (allocated(distribfft_arg%tab_fftdp2_local)) then
455   ABI_DEALLOCATE(distribfft_arg%tab_fftdp2_local)
456  end if
457  if (allocated(distribfft_arg%tab_fftdp3_local)) then
458   ABI_DEALLOCATE(distribfft_arg%tab_fftdp3_local)
459  end if
460  if (allocated(distribfft_arg%tab_fftwf2dg_local)) then
461   ABI_DEALLOCATE(distribfft_arg%tab_fftwf2dg_local)
462  end if
463  if (allocated(distribfft_arg%tab_fftdp2dg_local)) then
464    ABI_DEALLOCATE(distribfft_arg%tab_fftdp2dg_local)
465  end if
466  if (allocated(distribfft_arg%tab_fftdp3dg_local)) then
467    ABI_DEALLOCATE(distribfft_arg%tab_fftdp3dg_local)
468  end if
469 
470  DBG_EXIT("COLL")
471 
472 end subroutine destroy_distribfft

m_distribfft/distribfft_type [ Types ]

[ Top ] [ m_distribfft ] [ Types ]

NAME

 distribfft_type

FUNCTION

 The distribfft_type structured datatype gather different information
 for plane sharing for FFT parallelization

TODO

   1) One should create two separated tables: one for the wavefunctions and the other one for
      fourdp on the dense/coarse mesh.

   2) Use shorter names --> fftabs_type

SOURCE

 55  type, public :: distribfft_type
 56 
 57   integer :: n2_coarse=0
 58   ! Number of points along the y directions for the coarse FFT mesh
 59 
 60   integer :: n2_fine=0
 61   ! Number of points along the y directions for the dense FFT mesh
 62 
 63   !integer :: me_g0
 64   ! 1 if this MPI node has G=0, 0 otherwise.
 65   ! Needed for the FFTs of the wavefunctions.
 66 
 67   integer, allocatable :: tab_fftwf2_distrib(:)
 68   ! rank of the processors which own fft planes in 2nd dimension for fourwf
 69 
 70   integer, allocatable :: tab_fftdp2_distrib(:)
 71   ! rank of the processors which own fft planes in 2nd dimension for fourdp
 72 
 73   integer, allocatable :: tab_fftdp3_distrib(:)
 74   ! rank of the processors which own fft planes in 3rd dimension for fourdp
 75 
 76   integer, allocatable :: tab_fftwf2dg_distrib(:)
 77   ! rank of the processors which own fft planes in 2nd dimension for fourwf on fine grid
 78 
 79   integer, allocatable :: tab_fftdp2dg_distrib(:)
 80   ! rank of the processors which own fft planes in 2nd dimension for fourdp on fine grid
 81 
 82   integer, allocatable :: tab_fftdp3dg_distrib(:)
 83   ! rank of the processors which own fft planes in 3rd dimension for fourdp on fine grid
 84 
 85   integer, allocatable :: tab_fftwf2_local(:)
 86   ! local i2 indices in fourwf
 87 
 88   integer, allocatable :: tab_fftdp2_local(:)
 89   ! local i2 indices in fourdp
 90 
 91   integer, allocatable :: tab_fftdp3_local(:)
 92   ! local i3 indices in fourdp
 93 
 94   integer, allocatable :: tab_fftwf2dg_local(:)
 95   ! local i2 indices in fourwf on fine grid
 96 
 97   integer, allocatable :: tab_fftdp2dg_local(:)
 98   ! local i2 indices in fourdp on fine grid
 99 
100   integer, allocatable :: tab_fftdp3dg_local(:)
101   ! local i3 indices in fourdp on fine grid
102 
103 end type distribfft_type
104 
105  public :: init_distribfft         ! Initializes mpi information for FFT distribution.
106  public :: init_distribfft_seq     ! Initializes a sequential FFT distribution.
107  public :: destroy_distribfft      ! Free dynamic memory.
108  public :: copy_distribfft         ! Copy datatype.

m_distribfft/init_distribfft [ Functions ]

[ Top ] [ m_distribfft ] [ Functions ]

NAME

  init_distribfft

FUNCTION

  Initializes mpi information for FFT distribution
  Note that we always use cyclic distribution mode for the wavefunctions in G-space.
  MPI-FFT routines should always be compatible with this distribution.

INPUTS

 grid_type = 'c' or 'f' for information about coarse or fine fft grid
 nproc_fft = number of process used to distribute the fft
 n2,n3     = sizes of second and third fft grid

SIDE EFFECTS

  distribfft = instance of distribfft_type to initialize
  Update of "fft distrib" tabs accordingly to the fft parallelisation

PARENTS

      m_fft,m_fft_prof,m_qparticles,m_wfd,mpi_setup,scfcv,vtorhorec

CHILDREN

SOURCE

139 subroutine init_distribfft(distribfft_arg,grid_type,nproc_fft,n2,n3)
140 
141 
142 !This section has been created automatically by the script Abilint (TD).
143 !Do not modify the following lines by hand.
144 #undef ABI_FUNC
145 #define ABI_FUNC 'init_distribfft'
146 !End of the abilint section
147 
148  implicit none
149 
150 !Arguments ------------------------------------
151 !scalars
152  integer, intent(in) :: nproc_fft,n2,n3
153  character(len=1),intent(in) :: grid_type
154  type(distribfft_type), intent(inout) :: distribfft_arg
155 
156 !Local variables-------------------------------
157 !scalars
158  integer :: i2,i3,n2_local,n3_local
159  !character(len=500) :: msg
160 
161 ! ***********************************************************************
162 
163  DBG_ENTER("COLL")
164 
165  !local sizes
166  n2_local = n2 / nproc_fft
167  n3_local = n3 / nproc_fft
168 
169  select case (grid_type)
170  case ('c')
171     ! Updating information about coarse fft grid
172     if(distribfft_arg%n2_coarse > 0) then
173       if(n2 == distribfft_arg%n2_coarse) then
174         MSG_WARNING("The distribfft passed was already allocated for coarse grid on the same size")
175         return
176       else
177         MSG_ERROR("The distribfft passed was already allocated for coarse grid")
178       endif
179     end if
180     distribfft_arg%n2_coarse = n2
181     ! Initialisation of fft distrib tab
182     ABI_ALLOCATE(distribfft_arg%tab_fftwf2_distrib,(n2))
183     ABI_ALLOCATE(distribfft_arg%tab_fftwf2_local,(n2))
184     ABI_ALLOCATE(distribfft_arg%tab_fftdp2_distrib,(n2))
185     ABI_ALLOCATE(distribfft_arg%tab_fftdp2_local,(n2))
186     ABI_ALLOCATE(distribfft_arg%tab_fftdp3_distrib,(n3))
187     ABI_ALLOCATE(distribfft_arg%tab_fftdp3_local,(n3))
188     do i2=1, n2
189       ! Cyclic distribution of ig2 planes over fft processors
190       distribfft_arg%tab_fftwf2_distrib(i2) = modulo((i2-1),nproc_fft)
191       distribfft_arg%tab_fftwf2_local(i2)    = (i2-1)/nproc_fft + 1
192       ! Block distribution of i2 planes over fft processors for fourdp
193       distribfft_arg%tab_fftdp2_distrib(i2) = (i2-1) /  n2_local
194       distribfft_arg%tab_fftdp2_local(i2)   = modulo((i2-1),n2_local) + 1
195     end do
196     do i3=1, n3
197       ! Block distribution of i3 planes over fft processors for fourdp
198       distribfft_arg%tab_fftdp3_distrib(i3) = (i3-1) /  n3_local
199       distribfft_arg%tab_fftdp3_local(i3)   = modulo((i3-1),n3_local) + 1
200     end do
201 
202  case ('f')
203     if(distribfft_arg%n2_fine > 0) then
204       if(n2 == distribfft_arg%n2_fine) then
205         MSG_WARNING("The distribfft passed was already allocated for fine grid on the same size")
206         return
207       else
208         MSG_ERROR("The distribfft passed was already allocated for fine grid")
209       end if
210     endif
211     distribfft_arg%n2_fine = n2
212     ! Updating information about fine fft grid
213     ABI_ALLOCATE(distribfft_arg%tab_fftwf2dg_distrib,(n2))
214     ABI_ALLOCATE(distribfft_arg%tab_fftwf2dg_local,(n2))
215     ABI_ALLOCATE(distribfft_arg%tab_fftdp2dg_distrib,(n2))
216     ABI_ALLOCATE(distribfft_arg%tab_fftdp2dg_local,(n2))
217     ABI_ALLOCATE(distribfft_arg%tab_fftdp3dg_distrib,(n3))
218     ABI_ALLOCATE(distribfft_arg%tab_fftdp3dg_local,(n3))
219     do i2=1, n2
220       ! Cyclic distribution of ig2 planes over fft processors
221       distribfft_arg%tab_fftwf2dg_distrib(i2) = modulo((i2-1),nproc_fft)
222       distribfft_arg%tab_fftwf2dg_local(i2)    = (i2-1)/nproc_fft + 1
223       ! Block distribution of i2 planes over fft processors for fourdp on fine grid
224       distribfft_arg%tab_fftdp2dg_distrib(i2) = (i2-1) /  n2_local
225       distribfft_arg%tab_fftdp2dg_local(i2)   = modulo((i2-1),n2_local) + 1
226     end do
227     do i3=1, n3
228       ! Block distribution of i3 planes over fft processors for fourdp on fine grid
229       distribfft_arg%tab_fftdp3dg_distrib(i3) = (i3-1) /  n3_local
230       distribfft_arg%tab_fftdp3dg_local(i3)   = modulo((i3-1),n3_local) + 1
231     end do
232 
233  case default
234     MSG_ERROR("Unknown kind of fft grid! Only 'c' for coarse grid and 'f' for fine grid are allowed")
235  end select
236 
237  ! One needs to know if this node has G=0 when we do the FFTs of the wavefunctions
238  !distribfft_arg%me_g0 = 0
239  !if (distribfft_arg%tab_fftwf2_distrib(1) == me_fft) distribfft_arg%me_g0 = 1
240 
241  DBG_EXIT("COLL")
242 
243 end subroutine init_distribfft

m_distribfft/init_distribfft_seq [ Functions ]

[ Top ] [ m_distribfft ] [ Functions ]

NAME

  init_distribfft_seq

FUNCTION

  Initializes a sequential FFT distribution

INPUTS

 grid_type  = 'c' or 'f' for information about coarse or fine fft grid
 n2,n3 = sizes of second and third fft grid
 type_four = 'fourdp' or 'fourwf' or 'all' to prepare a call to fourdp/fourwf

SIDE EFFECTS

  distribfft = instance of t_distribfft to initialize
  Update of "fft distrib" tabs accordingly to the fft parallelisation

PARENTS

      atm2fft,bethe_salpeter,calc_vhxc_me,cut3d,dfpt_atm2fft,dieltcel,eph
      ks_ddiago,m_cut3d,m_dvdb,m_fft_prof,m_gsphere,m_ioarr,m_kxc,m_ppmodel
      m_screening,m_wfk,multipoles_fftr,pawgrnl,pawmknhat,pawmknhat_psipsi
      pawsushat,scfcv,screening,sigma,suscep_stat,susk,suskmm,wfk_analyze

CHILDREN

SOURCE

274 subroutine init_distribfft_seq(distribfft_arg,grid_type,n2,n3,type_four)
275 
276 
277 !This section has been created automatically by the script Abilint (TD).
278 !Do not modify the following lines by hand.
279 #undef ABI_FUNC
280 #define ABI_FUNC 'init_distribfft_seq'
281 !End of the abilint section
282 
283  implicit none
284 
285 
286 !Arguments ------------------------------------
287 !scalars
288  integer, intent(in) :: n2,n3
289  character(len=1),intent(in) :: grid_type
290  character(len=*),intent(in) :: type_four
291  type(distribfft_type), intent(inout) :: distribfft_arg
292 
293 !Local variables-------------------------------
294 !scalars
295  integer :: ii
296 
297 ! ***********************************************************************
298 
299  DBG_ENTER("COLL")
300 
301  !distribfft_arg%me_g0 = 1
302 
303  select case (grid_type)
304  case ('c')
305    distribfft_arg%n2_coarse = n2
306    if (type_four=='fourwf'.or.type_four(1:3)=='all') then
307      if (allocated(distribfft_arg%tab_fftwf2_distrib)) then
308        ABI_DEALLOCATE(distribfft_arg%tab_fftwf2_distrib)
309      end if
310      if (allocated(distribfft_arg%tab_fftwf2_local)) then
311        ABI_DEALLOCATE(distribfft_arg%tab_fftwf2_local)
312      end if
313      ABI_ALLOCATE(distribfft_arg%tab_fftwf2_distrib,(n2))
314      ABI_ALLOCATE(distribfft_arg%tab_fftwf2_local,(n2))
315      distribfft_arg%tab_fftwf2_distrib=0
316      distribfft_arg%tab_fftwf2_local=(/(ii,ii=1,n2)/)
317    end if
318    if (type_four=='fourdp'.or.type_four(1:3)=='all') then
319      if (allocated(distribfft_arg%tab_fftdp2_distrib)) then
320        ABI_DEALLOCATE(distribfft_arg%tab_fftdp2_distrib)
321      end if
322      if (allocated(distribfft_arg%tab_fftdp2_local)) then
323        ABI_DEALLOCATE(distribfft_arg%tab_fftdp2_local)
324      end if
325      if (allocated(distribfft_arg%tab_fftdp3_distrib)) then
326        ABI_DEALLOCATE(distribfft_arg%tab_fftdp3_distrib)
327      end if
328      if (allocated(distribfft_arg%tab_fftdp3_local)) then
329        ABI_DEALLOCATE(distribfft_arg%tab_fftdp3_local)
330      end if
331      ABI_ALLOCATE(distribfft_arg%tab_fftdp2_distrib,(n2))
332      ABI_ALLOCATE(distribfft_arg%tab_fftdp2_local,(n2))
333      ABI_ALLOCATE(distribfft_arg%tab_fftdp3_distrib,(n3))
334      ABI_ALLOCATE(distribfft_arg%tab_fftdp3_local,(n3))
335      distribfft_arg%tab_fftdp2_distrib=0
336      distribfft_arg%tab_fftdp3_distrib=0
337      distribfft_arg%tab_fftdp2_local=(/(ii,ii=1,n2)/)
338      distribfft_arg%tab_fftdp3_local=(/(ii,ii=1,n3)/)
339    end if
340 
341  case ('f')
342    distribfft_arg%n2_fine = n2
343    if (type_four=='fourwf'.or.type_four(1:3)=='all') then
344      if (allocated(distribfft_arg%tab_fftwf2dg_distrib)) then
345        ABI_DEALLOCATE(distribfft_arg%tab_fftwf2dg_distrib)
346      end if
347      if (allocated(distribfft_arg%tab_fftwf2dg_local)) then
348        ABI_DEALLOCATE(distribfft_arg%tab_fftwf2dg_local)
349      end if
350      ABI_ALLOCATE(distribfft_arg%tab_fftwf2dg_distrib,(n2))
351      ABI_ALLOCATE(distribfft_arg%tab_fftwf2dg_local,(n2))
352      distribfft_arg%tab_fftwf2dg_distrib=0
353      distribfft_arg%tab_fftwf2dg_local=(/(ii,ii=1,n2)/)
354    end if
355    if (type_four=='fourdp'.or.type_four(1:3)=='all') then
356      if (allocated(distribfft_arg%tab_fftdp2dg_distrib)) then
357        ABI_DEALLOCATE(distribfft_arg%tab_fftdp2dg_distrib)
358      end if
359      if (allocated(distribfft_arg%tab_fftdp2dg_local)) then
360        ABI_DEALLOCATE(distribfft_arg%tab_fftdp2dg_local)
361      end if
362      if (allocated(distribfft_arg%tab_fftdp3dg_distrib)) then
363        ABI_DEALLOCATE(distribfft_arg%tab_fftdp3dg_distrib)
364      end if
365      if (allocated(distribfft_arg%tab_fftdp3dg_local)) then
366        ABI_DEALLOCATE(distribfft_arg%tab_fftdp3dg_local)
367      end if
368      ABI_ALLOCATE(distribfft_arg%tab_fftdp2dg_distrib,(n2))
369      ABI_ALLOCATE(distribfft_arg%tab_fftdp2dg_local,(n2))
370      ABI_ALLOCATE(distribfft_arg%tab_fftdp3dg_distrib,(n3))
371      ABI_ALLOCATE(distribfft_arg%tab_fftdp3dg_local,(n3))
372      distribfft_arg%tab_fftdp2dg_distrib=0
373      distribfft_arg%tab_fftdp3dg_distrib=0
374      distribfft_arg%tab_fftdp2dg_local=(/(ii,ii=1,n2)/)
375      distribfft_arg%tab_fftdp3dg_local=(/(ii,ii=1,n3)/)
376    end if
377 
378  case default
379     MSG_ERROR("Unknown kind of fft grid! Only 'c' for coarse grid and 'f' for fine grid are allowed")
380  end select
381 
382  DBG_EXIT("COLL")
383 
384 end subroutine init_distribfft_seq