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-2022 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 .

SOURCE

17 #if defined HAVE_CONFIG_H
18 #include "config.h"
19 #endif
20 
21 #include "abi_common.h"
22 
23 MODULE m_distribfft
24 
25  use defs_basis
26  use m_abicore
27  use m_errors
28 
29  implicit none
30 
31  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

SOURCE

437 subroutine copy_distribfft(distribfft_src, distribfft_dst)
438 
439 !Arguments ------------------------------------
440  type(distribfft_type),intent(in)   :: distribfft_src
441  type(distribfft_type),intent(out) :: distribfft_dst
442 
443 ! ***********************************************************************
444 
445  DBG_ENTER("COLL")
446 
447  distribfft_dst%n2_coarse=  distribfft_src%n2_coarse
448  distribfft_dst%n2_fine  =  distribfft_src%n2_fine
449 
450  if (allocated(distribfft_src%tab_fftwf2_distrib)) then
451    ABI_MALLOC(distribfft_dst%tab_fftwf2_distrib,(size(distribfft_src%tab_fftwf2_distrib)))
452    distribfft_dst%tab_fftwf2_distrib=distribfft_src%tab_fftwf2_distrib
453  end if
454 
455  if (allocated(distribfft_src%tab_fftdp2_distrib)) then
456   ABI_MALLOC(distribfft_dst%tab_fftdp2_distrib,(size(distribfft_src%tab_fftdp2_distrib)))
457   distribfft_dst%tab_fftdp2_distrib=distribfft_src%tab_fftdp2_distrib
458  end if
459 
460  if (allocated(distribfft_src%tab_fftdp3_distrib)) then
461   ABI_MALLOC(distribfft_dst%tab_fftdp3_distrib,(size(distribfft_src%tab_fftdp3_distrib)))
462   distribfft_dst%tab_fftdp3_distrib=distribfft_src%tab_fftdp3_distrib
463  end if
464 
465  if (allocated(distribfft_src%tab_fftwf2dg_distrib)) then
466    ABI_MALLOC(distribfft_dst%tab_fftwf2dg_distrib,(size(distribfft_src%tab_fftwf2dg_distrib)))
467    distribfft_dst%tab_fftwf2dg_distrib=distribfft_src%tab_fftwf2dg_distrib
468  end if
469 
470  if (allocated(distribfft_src%tab_fftdp2dg_distrib)) then
471    ABI_MALLOC(distribfft_dst%tab_fftdp2dg_distrib,(size(distribfft_src%tab_fftdp2dg_distrib)))
472    distribfft_dst%tab_fftdp2dg_distrib=distribfft_src%tab_fftdp2dg_distrib
473  end if
474 
475  if (allocated(distribfft_src%tab_fftdp3dg_distrib)) then
476    ABI_MALLOC(distribfft_dst%tab_fftdp3dg_distrib,(size(distribfft_src%tab_fftdp3dg_distrib)))
477    distribfft_dst%tab_fftdp3dg_distrib=distribfft_src%tab_fftdp3dg_distrib
478  end if
479 
480  if (allocated(distribfft_src%tab_fftwf2_local)) then
481    ABI_MALLOC(distribfft_dst%tab_fftwf2_local,(size(distribfft_src%tab_fftwf2_local)))
482    distribfft_dst%tab_fftwf2_local=distribfft_src%tab_fftwf2_local
483  end if
484 
485  if (allocated(distribfft_src%tab_fftdp2_local)) then
486    ABI_MALLOC(distribfft_dst%tab_fftdp2_local,(size(distribfft_src%tab_fftdp2_local)))
487    distribfft_dst%tab_fftdp2_local=distribfft_src%tab_fftdp2_local
488  end if
489 
490  if (allocated(distribfft_src%tab_fftdp3_local)) then
491    ABI_MALLOC(distribfft_dst%tab_fftdp3_local,(size(distribfft_src%tab_fftdp3_local)))
492    distribfft_dst%tab_fftdp3_local=distribfft_src%tab_fftdp3_local
493  end if
494 
495  if (allocated(distribfft_src%tab_fftwf2dg_local)) then
496    ABI_MALLOC(distribfft_dst%tab_fftwf2dg_local,(size(distribfft_src%tab_fftwf2dg_local)))
497    distribfft_dst%tab_fftwf2dg_local=distribfft_src%tab_fftwf2dg_local
498  end if
499 
500  if (allocated(distribfft_src%tab_fftdp2dg_local)) then
501    ABI_MALLOC(distribfft_dst%tab_fftdp2dg_local,(size(distribfft_src%tab_fftdp2dg_local)))
502    distribfft_dst%tab_fftdp2dg_local=distribfft_src%tab_fftdp2dg_local
503  end if
504 
505  if (allocated(distribfft_src%tab_fftdp3dg_local)) then
506    ABI_MALLOC(distribfft_dst%tab_fftdp3dg_local,(size(distribfft_src%tab_fftdp3dg_local)))
507    distribfft_dst%tab_fftdp3dg_local=distribfft_src%tab_fftdp3dg_local
508  end if
509 
510  DBG_EXIT("COLL")
511 
512 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

SOURCE

365 subroutine destroy_distribfft(distribfft_arg)
366 
367 !Arguments ------------------------------------
368  type(distribfft_type), intent(inout) :: distribfft_arg
369 
370 ! ***********************************************************************
371 
372  DBG_ENTER("COLL")
373 
374  distribfft_arg%n2_coarse=0
375  distribfft_arg%n2_fine  =0
376 
377  if (allocated(distribfft_arg%tab_fftwf2_distrib)) then
378    ABI_FREE(distribfft_arg%tab_fftwf2_distrib)
379  end if
380 
381  if (allocated(distribfft_arg%tab_fftdp2_distrib)) then
382    ABI_FREE(distribfft_arg%tab_fftdp2_distrib)
383  end if
384  if (allocated(distribfft_arg%tab_fftdp3_distrib)) then
385    ABI_FREE(distribfft_arg%tab_fftdp3_distrib)
386  end if
387  if (allocated(distribfft_arg%tab_fftwf2dg_distrib)) then
388   ABI_FREE(distribfft_arg%tab_fftwf2dg_distrib)
389  end if
390  if (allocated(distribfft_arg%tab_fftdp2dg_distrib)) then
391   ABI_FREE(distribfft_arg%tab_fftdp2dg_distrib)
392  end if
393  if (allocated(distribfft_arg%tab_fftdp3dg_distrib)) then
394   ABI_FREE(distribfft_arg%tab_fftdp3dg_distrib)
395  end if
396 
397  if (allocated(distribfft_arg%tab_fftwf2_local)) then
398   ABI_FREE(distribfft_arg%tab_fftwf2_local)
399  end if
400 
401  if (allocated(distribfft_arg%tab_fftdp2_local)) then
402   ABI_FREE(distribfft_arg%tab_fftdp2_local)
403  end if
404  if (allocated(distribfft_arg%tab_fftdp3_local)) then
405   ABI_FREE(distribfft_arg%tab_fftdp3_local)
406  end if
407  if (allocated(distribfft_arg%tab_fftwf2dg_local)) then
408   ABI_FREE(distribfft_arg%tab_fftwf2dg_local)
409  end if
410  if (allocated(distribfft_arg%tab_fftdp2dg_local)) then
411    ABI_FREE(distribfft_arg%tab_fftdp2dg_local)
412  end if
413  if (allocated(distribfft_arg%tab_fftdp3dg_local)) then
414    ABI_FREE(distribfft_arg%tab_fftdp3dg_local)
415  end if
416 
417  DBG_EXIT("COLL")
418 
419 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

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

SOURCE

129 subroutine init_distribfft(distribfft_arg,grid_type,nproc_fft,n2,n3)
130 
131 !Arguments ------------------------------------
132 !scalars
133  integer, intent(in) :: nproc_fft,n2,n3
134  character(len=1),intent(in) :: grid_type
135  type(distribfft_type), intent(inout) :: distribfft_arg
136 
137 !Local variables-------------------------------
138 !scalars
139  integer :: i2,i3,n2_local,n3_local
140  !character(len=500) :: msg
141 
142 ! ***********************************************************************
143 
144  DBG_ENTER("COLL")
145 
146  !local sizes
147  n2_local = n2 / nproc_fft
148  n3_local = n3 / nproc_fft
149 
150  select case (grid_type)
151  case ('c')
152     ! Updating information about coarse fft grid
153     if(distribfft_arg%n2_coarse > 0) then
154       if(n2 == distribfft_arg%n2_coarse) then
155         ABI_WARNING("The distribfft passed was already allocated for coarse grid on the same size")
156         return
157       else
158         ABI_ERROR("The distribfft passed was already allocated for coarse grid")
159       endif
160     end if
161     distribfft_arg%n2_coarse = n2
162     ! Initialisation of fft distrib tab
163     ABI_MALLOC(distribfft_arg%tab_fftwf2_distrib,(n2))
164     ABI_MALLOC(distribfft_arg%tab_fftwf2_local,(n2))
165     ABI_MALLOC(distribfft_arg%tab_fftdp2_distrib,(n2))
166     ABI_MALLOC(distribfft_arg%tab_fftdp2_local,(n2))
167     ABI_MALLOC(distribfft_arg%tab_fftdp3_distrib,(n3))
168     ABI_MALLOC(distribfft_arg%tab_fftdp3_local,(n3))
169     do i2=1, n2
170       ! Cyclic distribution of ig2 planes over fft processors
171       distribfft_arg%tab_fftwf2_distrib(i2) = modulo((i2-1),nproc_fft)
172       distribfft_arg%tab_fftwf2_local(i2)    = (i2-1)/nproc_fft + 1
173       ! Block distribution of i2 planes over fft processors for fourdp
174       distribfft_arg%tab_fftdp2_distrib(i2) = (i2-1) /  n2_local
175       distribfft_arg%tab_fftdp2_local(i2)   = modulo((i2-1),n2_local) + 1
176     end do
177     do i3=1, n3
178       ! Block distribution of i3 planes over fft processors for fourdp
179       distribfft_arg%tab_fftdp3_distrib(i3) = (i3-1) /  n3_local
180       distribfft_arg%tab_fftdp3_local(i3)   = modulo((i3-1),n3_local) + 1
181     end do
182 
183  case ('f')
184     if(distribfft_arg%n2_fine > 0) then
185       if(n2 == distribfft_arg%n2_fine) then
186         ABI_WARNING("The distribfft passed was already allocated for fine grid on the same size")
187         return
188       else
189         ABI_ERROR("The distribfft passed was already allocated for fine grid")
190       end if
191     endif
192     distribfft_arg%n2_fine = n2
193     ! Updating information about fine fft grid
194     ABI_MALLOC(distribfft_arg%tab_fftwf2dg_distrib,(n2))
195     ABI_MALLOC(distribfft_arg%tab_fftwf2dg_local,(n2))
196     ABI_MALLOC(distribfft_arg%tab_fftdp2dg_distrib,(n2))
197     ABI_MALLOC(distribfft_arg%tab_fftdp2dg_local,(n2))
198     ABI_MALLOC(distribfft_arg%tab_fftdp3dg_distrib,(n3))
199     ABI_MALLOC(distribfft_arg%tab_fftdp3dg_local,(n3))
200     do i2=1, n2
201       ! Cyclic distribution of ig2 planes over fft processors
202       distribfft_arg%tab_fftwf2dg_distrib(i2) = modulo((i2-1),nproc_fft)
203       distribfft_arg%tab_fftwf2dg_local(i2)    = (i2-1)/nproc_fft + 1
204       ! Block distribution of i2 planes over fft processors for fourdp on fine grid
205       distribfft_arg%tab_fftdp2dg_distrib(i2) = (i2-1) /  n2_local
206       distribfft_arg%tab_fftdp2dg_local(i2)   = modulo((i2-1),n2_local) + 1
207     end do
208     do i3=1, n3
209       ! Block distribution of i3 planes over fft processors for fourdp on fine grid
210       distribfft_arg%tab_fftdp3dg_distrib(i3) = (i3-1) /  n3_local
211       distribfft_arg%tab_fftdp3dg_local(i3)   = modulo((i3-1),n3_local) + 1
212     end do
213 
214  case default
215     ABI_ERROR("Unknown kind of fft grid! Only 'c' for coarse grid and 'f' for fine grid are allowed")
216  end select
217 
218  ! One needs to know if this node has G=0 when we do the FFTs of the wavefunctions
219  !distribfft_arg%me_g0 = 0
220  !if (distribfft_arg%tab_fftwf2_distrib(1) == me_fft) distribfft_arg%me_g0 = 1
221 
222  DBG_EXIT("COLL")
223 
224 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

SOURCE

247 subroutine init_distribfft_seq(distribfft_arg,grid_type,n2,n3,type_four)
248 
249 !Arguments ------------------------------------
250 !scalars
251  integer, intent(in) :: n2,n3
252  character(len=1),intent(in) :: grid_type
253  character(len=*),intent(in) :: type_four
254  type(distribfft_type), intent(inout) :: distribfft_arg
255 
256 !Local variables-------------------------------
257 !scalars
258  integer :: ii
259 
260 ! ***********************************************************************
261 
262  DBG_ENTER("COLL")
263 
264  !distribfft_arg%me_g0 = 1
265 
266  select case (grid_type)
267  case ('c')
268    distribfft_arg%n2_coarse = n2
269    if (type_four=='fourwf'.or.type_four(1:3)=='all') then
270      if (allocated(distribfft_arg%tab_fftwf2_distrib)) then
271        ABI_FREE(distribfft_arg%tab_fftwf2_distrib)
272      end if
273      if (allocated(distribfft_arg%tab_fftwf2_local)) then
274        ABI_FREE(distribfft_arg%tab_fftwf2_local)
275      end if
276      ABI_MALLOC(distribfft_arg%tab_fftwf2_distrib,(n2))
277      ABI_MALLOC(distribfft_arg%tab_fftwf2_local,(n2))
278      distribfft_arg%tab_fftwf2_distrib=0
279      distribfft_arg%tab_fftwf2_local=(/(ii,ii=1,n2)/)
280    end if
281    if (type_four=='fourdp'.or.type_four(1:3)=='all') then
282      if (allocated(distribfft_arg%tab_fftdp2_distrib)) then
283        ABI_FREE(distribfft_arg%tab_fftdp2_distrib)
284      end if
285      if (allocated(distribfft_arg%tab_fftdp2_local)) then
286        ABI_FREE(distribfft_arg%tab_fftdp2_local)
287      end if
288      if (allocated(distribfft_arg%tab_fftdp3_distrib)) then
289        ABI_FREE(distribfft_arg%tab_fftdp3_distrib)
290      end if
291      if (allocated(distribfft_arg%tab_fftdp3_local)) then
292        ABI_FREE(distribfft_arg%tab_fftdp3_local)
293      end if
294      ABI_MALLOC(distribfft_arg%tab_fftdp2_distrib,(n2))
295      ABI_MALLOC(distribfft_arg%tab_fftdp2_local,(n2))
296      ABI_MALLOC(distribfft_arg%tab_fftdp3_distrib,(n3))
297      ABI_MALLOC(distribfft_arg%tab_fftdp3_local,(n3))
298      distribfft_arg%tab_fftdp2_distrib=0
299      distribfft_arg%tab_fftdp3_distrib=0
300      distribfft_arg%tab_fftdp2_local=(/(ii,ii=1,n2)/)
301      distribfft_arg%tab_fftdp3_local=(/(ii,ii=1,n3)/)
302    end if
303 
304  case ('f')
305    distribfft_arg%n2_fine = n2
306    if (type_four=='fourwf'.or.type_four(1:3)=='all') then
307      if (allocated(distribfft_arg%tab_fftwf2dg_distrib)) then
308        ABI_FREE(distribfft_arg%tab_fftwf2dg_distrib)
309      end if
310      if (allocated(distribfft_arg%tab_fftwf2dg_local)) then
311        ABI_FREE(distribfft_arg%tab_fftwf2dg_local)
312      end if
313      ABI_MALLOC(distribfft_arg%tab_fftwf2dg_distrib,(n2))
314      ABI_MALLOC(distribfft_arg%tab_fftwf2dg_local,(n2))
315      distribfft_arg%tab_fftwf2dg_distrib=0
316      distribfft_arg%tab_fftwf2dg_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_fftdp2dg_distrib)) then
320        ABI_FREE(distribfft_arg%tab_fftdp2dg_distrib)
321      end if
322      if (allocated(distribfft_arg%tab_fftdp2dg_local)) then
323        ABI_FREE(distribfft_arg%tab_fftdp2dg_local)
324      end if
325      if (allocated(distribfft_arg%tab_fftdp3dg_distrib)) then
326        ABI_FREE(distribfft_arg%tab_fftdp3dg_distrib)
327      end if
328      if (allocated(distribfft_arg%tab_fftdp3dg_local)) then
329        ABI_FREE(distribfft_arg%tab_fftdp3dg_local)
330      end if
331      ABI_MALLOC(distribfft_arg%tab_fftdp2dg_distrib,(n2))
332      ABI_MALLOC(distribfft_arg%tab_fftdp2dg_local,(n2))
333      ABI_MALLOC(distribfft_arg%tab_fftdp3dg_distrib,(n3))
334      ABI_MALLOC(distribfft_arg%tab_fftdp3dg_local,(n3))
335      distribfft_arg%tab_fftdp2dg_distrib=0
336      distribfft_arg%tab_fftdp3dg_distrib=0
337      distribfft_arg%tab_fftdp2dg_local=(/(ii,ii=1,n2)/)
338      distribfft_arg%tab_fftdp3dg_local=(/(ii,ii=1,n3)/)
339    end if
340 
341  case default
342     ABI_ERROR("Unknown kind of fft grid! Only 'c' for coarse grid and 'f' for fine grid are allowed")
343  end select
344 
345  DBG_EXIT("COLL")
346 
347 end subroutine init_distribfft_seq