MONC
Data Types | Functions/Subroutines | Variables
pencil_fft_mod Module Reference

This Pencil FFT performs 3D forward and backwards FFTs using pencil decomposition. It uses FFTW for the actual FFT kernel and this module contains all the data decomposition around this. There is no FFT required in Z, so this performs FFTs in Y and X (in that order forward and reversed backwards.) The data decomposition is the complex aspect, there is the concept of forward and backwards transformations. Forward transformations will go from pencil Z to Y to X and the backwards transformations undo these, so go from X to Y to Z. Note that we use quite a lot of buffer space here, this could be cut down if Y=X dimensions so some optimisation on memory could be done there in that case. More...

Data Types

type  pencil_transposition
 Describes a specific pencil transposition, from one pencil decomposition to another. More...
 

Functions/Subroutines

integer function, dimension(3), public initialise_pencil_fft (current_state, my_y_start, my_x_start)
 Initialises the pencil FFT functionality, this will create the transposition structures needed. More...
 
subroutine, public finalise_pencil_fft (monc_communicator)
 Cleans up allocated buffer memory. More...
 
subroutine, public perform_forward_3dfft (current_state, source_data, target_data)
 Performs a forward 3D FFT and currently results in target data which is the X, Z, Y oriented pencil Note that the source_data here takes no account for the halo, it is up to caller to exclude this. This does no FFT in Z, but transposes to Y, does FFT in Y, then transposes to X and performs an FFT in that dimension. Pencil decomposition is used which has already been set up. More...
 
subroutine, public perform_backwards_3dfft (current_state, source_data, target_data)
 Performs a backwards 3D FFT and currently results in target data which is the X, Z, Y oriented pencil Note that the source_data here takes no account for the halo, it is up to caller to exclude this. This does no FFT in Z, but transposes to Y, does FFT in Y, then transposes to X and performs an FFT in that dimension. Pencil decomposition is used which has already been set up. More...
 
subroutine initialise_buffers ()
 Initialises memory for the buffers used in the FFT. More...
 
subroutine initialise_transpositions (current_state, y_distinct_sizes, x_distinct_sizes)
 Initialises the pencil transpositions, from a pencil in one dimension to that in another. More...
 
type(pencil_transposition) function create_transposition (global_grid, existing_transposition, new_pencil_dim, process_dim_sizes, direction, extended_dimensions)
 Creates a specific pencil transposition description. It is maybe more a decomposition description, but the main complexity comes from the transposition from existing decomposition to new decomposition so therefore it is called transposition. The new pencil decomposition depends not only on the dimension to split on, but also the existing pencil decomposition. The new decomposed dimension (i.e. the existing pencil dimension) * other local dimensions is used as the sending size, receiving though requires knowledge about the data size on the source process so others will send their pencil dimension size to this process. More...
 
subroutine transpose_and_forward_fft_in_y (current_state, source_data, buffer, real_buffer)
 Performs the transposition and forward FFT in the y dimension then converts back to real numbers. The Y size is (n/2+1)*2 due to the complex to real transformation after the FFT. More...
 
subroutine transpose_and_backward_fft_in_x (current_state, source_data, buffer, real_buffer)
 Performs the backwards FFT in X and then transposes to Y pencil. The FFT requires complex numbers which are converted to real, so the this real to complex operation is performed first. If n is the logical size of the FFT row, then the input size is n+2, complex number size is n/2+1 and we get n reals out. More...
 
subroutine transpose_and_forward_fft_in_x (current_state, buffer1, buffer2, buffer3)
 Performs the transposition and forward FFT in the x dimension. After the FFT the complex space is converted back into real numbers. The X size is (n/2+1)*2 due to this transformation. More...
 
subroutine transpose_and_backward_fft_in_y (current_state, source_data, buffer, real_buffer)
 Performs the backwards FFT in Y and then transposes to Z pencil. The FFT requires complex numbers which are converted to real, so the this real to complex operation is performed first. If n is the logical size of the FFT row, then the input size is n+2, complex number size is n/2+1 and we get n reals out. More...
 
subroutine transpose_to_pencil (transposition_description, source_dims, communicator, direction, source_data, target_data)
 Transposes globally to a new pencil decomposition. This goes from the source dimensions a,b,c to b,c,a (forwards) or c,a,b (backwards.) It requires multiple steps, first the local data is transposed to c,b,a regardless of direction. then it is communicated via alltoall, each process then assembles its own b,c,a or c,a,b data via contiguising across blocks as the data layout is nonlinear. More...
 
subroutine contiguise_data (transposition_description, source_dims, direction, source_real_buffer, target_real_buffer)
 Contiguises from c,b,a to b,c,a (forwards) or c,a,b (backwards) where these are defined by the source_dims argument. It is not as simple as just swapping the required dimensions, as this is after the mpi alltoall and each block lies after the previous block running sequentially in a. More...
 
subroutine perform_r2c_fft (source_data, transformed_data, row_size, num_rows, plan_id)
 Actually performs a forward real to complex FFT. More...
 
subroutine perform_c2r_fft (source_data, transformed_data, row_size, num_rows, plan_id)
 Performs the complex to real (backwards) FFT. More...
 
subroutine rearrange_data_for_sending (real_source, real_target)
 Rearranges data for sending, transposing a,b,c into c,b,a . This is done as alltoall splits on dimension c so to go from one pencil to another we assume here that a is the existing pencil as it is contiguous. More...
 
subroutine determine_my_process_sizes_per_dim (existing_pencil_dim, existing_pencil_size, new_pencil_procs_per_dim, global_grid, extended_dimensions, specific_sizes_per_dim)
 Determines the number of elements to on my process per dimension which either need to be sent to (forwards transformation) or received from (backwards) each target process (in the row or column) This depends on the existing pencil decomposition, as effectively we are breaking that contigulity and decomposing it into n blocks in that dimension now (provided by new_pencil_procs_per_dim) More...
 
subroutine determine_offsets_from_size (source_sizes, determined_offsets)
 Simple helper function to deduce send or receive offsets from the sizes. More...
 
integer function, dimension(3) determine_pencil_process_dimensions (new_pencil_dim, existing_pencil_dim, existing_pencil_procs)
 Determines the number of processes in each dimension for the target decomposition. This depends heavily on the existing decomposition, as we basically contiguise our pencil dimension and decompose the existing pencil dimension. The third dimension remains unchanged. More...
 
integer function, dimension(3) determine_my_pencil_location (new_pencil_dim, existing_pencil_dim, existing_locations)
 Determines my location for each dimension in the new pencil decomposition. I.e. which block I am operating on. More...
 
subroutine concatenate_dimension_sizes (dims, concatenated_dim_sizes)
 Concatenates sizes in multiple dimensions for each target process (in a row or column) into a product of that. This represents all the dimension sizes per process. More...
 
subroutine determine_matching_process_dimensions (new_pencil_dim, existing_pencil_dim, proc_sizes, my_pencil_size, pencil_processes_per_dim, specific_sizes_per_dim)
 Determines the sizes per dimension on the matching process either to receive from (forward transposition) or send to (backwards transposition) each source process. Not only does this depend on the my pencil sizes, but it also depends on the amount of data that the source process has to send over. More...
 
type(pencil_transposition) function create_initial_transposition_description (current_state)
 Creates an initial transposition representation of the Z pencil that MONC is normally decomposed in. This is then fed into the create transposition procedure which will generate transpositions to other pencils. More...
 
integer function, dimension(3) determine_pencil_size (new_pencil_dim, pencil_process_layout, my_pencil_location, existing_transposition, global_grid, extended_dimensions)
 Deduces the size of my (local) pencil based upon the new decomposition. This depends heavily on the current pencil decomposition, the new pencil dimension is the global size, the existing pencil dimension becomes decomposed based on the number of processes in that dimension. The third dimension remains unchanged. More...
 
logical function is_extended_dimension (dimension, extended_dimensions)
 Determines whether or not the specific dimension is in the list of extended dimensions. More...
 
integer function, dimension(size(process_dim_sizes)) normal_to_extended_process_dim_sizes (process_dim_sizes)
 Transforms real process dimension sizes into their real after FFT complex->real transformation. The way this works is that it goes from n to (n/2+1)*2 numbers which is distributed amongst the processes deterministically. More...
 
subroutine convert_complex_to_real (complex_data, real_data)
 Converts complex representation to its real data counterpart and is called after each forward FFT. After a r2c FFT, there are n/2+1 complex numbers - which means that there will be more real numbers in Fourier space than are provided into the forward FFT call (due to the extra +1). Note that the real size n will always be complex size * 2 This always unpacks the complex dimension in the first dimension. More...
 
subroutine convert_real_to_complex (real_data, complex_data)
 Converts reals into their complex representation, this is called for backwards FFTs as we need to feed in complex numbers to force FFTW to do a backwards. It is a relatively simple transformation, as n goes into n/2 complex numbers and as this is the result of the convert_complex_to_real procedure, n always divides evenly. This is always applied to the first dimension of the real data. More...
 
integer function deduce_my_global_start (current_state, dimension)
 Determines my global start coordinate in Fourier space. This is required for cos y and cos x calculation which is fed into the tridiagonal solver. After the forward FFTs, each process has ((n/2+1)/p+r) * 2 elements, where p is the number of processes and r is the uneven process remainder (1 or 0 depending on p). Therefore some processes will have t elements, and some t-2 elements to feed into the solver. More...
 

Variables

integer, parameter forward =1
 
integer, parameter backward =2
 Transposition directions. More...
 
integer dim_y_comm
 
integer dim_x_comm
 Communicators for each dimension. More...
 
type(pencil_transpositiony_from_z_transposition
 
type(pencil_transpositionx_from_y_transposition
 
type(pencil_transpositiony_from_x_transposition
 
type(pencil_transpositionz_from_y_transposition
 
type(pencil_transpositiony_from_z_2_transposition
 
type(pencil_transpositionx_from_y_2_transposition
 
type(pencil_transpositiony_from_x_2_transposition
 
type(pencil_transpositionz_from_y_2_transposition
 
real(kind=default_precision), dimension(:,:,:), pointer, contiguous real_buffer1
 
real(kind=default_precision), dimension(:,:,:), pointer, contiguous real_buffer2
 
real(kind=default_precision), dimension(:,:,:), pointer, contiguous real_buffer3
 
real(kind=default_precision), dimension(:,:,:), pointer, contiguous fft_in_y_buffer
 
real(kind=default_precision), dimension(:,:,:), pointer, contiguous fft_in_x_buffer
 
complex(c_double_complex), dimension(:,:,:), pointer, contiguous buffer1
 
complex(c_double_complex), dimension(:,:,:), pointer, contiguous buffer2
 
type(c_ptr), dimension(4) fftw_plan
 
logical, dimension(4) fftw_plan_initialised =.false.
 

Detailed Description

This Pencil FFT performs 3D forward and backwards FFTs using pencil decomposition. It uses FFTW for the actual FFT kernel and this module contains all the data decomposition around this. There is no FFT required in Z, so this performs FFTs in Y and X (in that order forward and reversed backwards.) The data decomposition is the complex aspect, there is the concept of forward and backwards transformations. Forward transformations will go from pencil Z to Y to X and the backwards transformations undo these, so go from X to Y to Z. Note that we use quite a lot of buffer space here, this could be cut down if Y=X dimensions so some optimisation on memory could be done there in that case.

Function/Subroutine Documentation

◆ concatenate_dimension_sizes()

subroutine pencil_fft_mod::concatenate_dimension_sizes ( integer, dimension(:,:), intent(in)  dims,
integer, dimension(:), intent(inout)  concatenated_dim_sizes 
)
private

Concatenates sizes in multiple dimensions for each target process (in a row or column) into a product of that. This represents all the dimension sizes per process.

Parameters
dimsThe sizes, per dimension and per process that we will fold into target process

Definition at line 574 of file pencilfft.F90.

574  integer, dimension(:,:), intent(in) :: dims
575  integer, dimension(:), intent(inout) :: concatenated_dim_sizes
576 
577  integer :: i
578 
579  do i=1,size(dims, 2)
580  concatenated_dim_sizes(i)=product(dims(:,i))
581  end do
Here is the caller graph for this function:

◆ contiguise_data()

subroutine pencil_fft_mod::contiguise_data ( type(pencil_transposition), intent(in)  transposition_description,
integer, dimension(3), intent(in)  source_dims,
integer, intent(in)  direction,
real(kind=default_precision), dimension(:), intent(in)  source_real_buffer,
real(kind=default_precision), dimension(:,:,:), intent(out)  target_real_buffer 
)
private

Contiguises from c,b,a to b,c,a (forwards) or c,a,b (backwards) where these are defined by the source_dims argument. It is not as simple as just swapping the required dimensions, as this is after the mpi alltoall and each block lies after the previous block running sequentially in a.

Parameters
transposition_descriptionTransposition descriptor
source_dimsRepresentation a,b,c of source data, will contiguise to b,a,c
directionWhether we wish to contiguise forwards or backwards
source_real_bufferSource real data to transform
target_real_bufferTarget real data which is the result of the operation

Definition at line 390 of file pencilfft.F90.

390  integer, intent(in) :: source_dims(3), direction
391  type(pencil_transposition), intent(in) :: transposition_description
392  real(kind=DEFAULT_PRECISION), dimension(:), intent(in) :: source_real_buffer
393  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(out) :: target_real_buffer
394 
395  integer :: number_blocks, i, j, k, n, index_prefix, index_prefix_dim, block_offset, source_index
396 
397  number_blocks=size(transposition_description%recv_sizes)
398  index_prefix=0
399  block_offset=0
400  index_prefix_dim=merge(2,1, direction == forward)
401  do i=1,number_blocks
402  if (i .ge. 2) then
403  index_prefix=index_prefix+transposition_description%recv_dims(source_dims(index_prefix_dim), i-1)
404  block_offset=block_offset+transposition_description%recv_sizes(i-1)
405  end if
406  !Transformation is either cba -> bca (forward) or cab (backwards)
407  do j=1, transposition_description%recv_dims(source_dims(3), i) ! a
408  do k=1, transposition_description%recv_dims(source_dims(1), i) ! c
409  do n=1, transposition_description%recv_dims(source_dims(2), i) ! b
410  source_index=block_offset+(j-1)* transposition_description%recv_dims(source_dims(1), i)* &
411  transposition_description%recv_dims(source_dims(2), i)+ (n-1)* &
412  transposition_description%recv_dims(source_dims(1), i)+k
413  if (direction == forward) then
414  target_real_buffer(index_prefix+n, k, j)=source_real_buffer(source_index) ! bca
415  else
416  target_real_buffer(index_prefix+k, j, n)=source_real_buffer(source_index) ! cab
417  end if
418  end do
419  end do
420  end do
421  end do
Here is the caller graph for this function:

◆ convert_complex_to_real()

subroutine pencil_fft_mod::convert_complex_to_real ( complex(c_double_complex), dimension(:,:,:), intent(in)  complex_data,
real(kind=default_precision), dimension(:,:,:), intent(out)  real_data 
)
private

Converts complex representation to its real data counterpart and is called after each forward FFT. After a r2c FFT, there are n/2+1 complex numbers - which means that there will be more real numbers in Fourier space than are provided into the forward FFT call (due to the extra +1). Note that the real size n will always be complex size * 2 This always unpacks the complex dimension in the first dimension.

Parameters
complex_dataComplex data in Z,Y,X orientation to be unpacked into its real representation
real_dataThe real representation is written into here

Definition at line 705 of file pencilfft.F90.

705  complex(C_DOUBLE_COMPLEX), dimension(:,:,:), intent(in) :: complex_data
706  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(out) :: real_data
707 
708  integer :: i, j, k
709 
710  do i=1,size(real_data,3)
711  do j=1,size(real_data,2)
712  do k=1,size(real_data,1),2
713  real_data(k,j,i)=real(real(complex_data((k+1)/2,j,i)), kind=default_precision)
714  real_data(k+1,j,i)=real(aimag(complex_data((k+1)/2,j,i)), kind=default_precision)
715  end do
716  end do
717  end do
Here is the caller graph for this function:

◆ convert_real_to_complex()

subroutine pencil_fft_mod::convert_real_to_complex ( real(kind=default_precision), dimension(:,:,:), intent(in)  real_data,
complex(c_double_complex), dimension(:,:,:), intent(out), pointer, contiguous  complex_data 
)
private

Converts reals into their complex representation, this is called for backwards FFTs as we need to feed in complex numbers to force FFTW to do a backwards. It is a relatively simple transformation, as n goes into n/2 complex numbers and as this is the result of the convert_complex_to_real procedure, n always divides evenly. This is always applied to the first dimension of the real data.

Parameters
real_dataThe source real data to pack into the complex data, it is oriented Z,Y,X
complex_dataTarget complex data which the real data is packaged into

Definition at line 727 of file pencilfft.F90.

727  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(in) :: real_data
728  complex(C_DOUBLE_COMPLEX), dimension(:,:,:), contiguous, pointer, intent(out) :: complex_data
729 
730  integer :: i, j, k
731 
732  complex_data=cmplx(0.0d0, 0.0d0, kind=c_double_complex)
733 
734  do i=1,size(real_data,3)
735  do j=1,size(real_data,2)
736  do k=1,size(real_data,1),2
737  complex_data((k+1)/2,j,i)=cmplx(real_data(k,j,i), real_data(k+1,j,i), kind=c_double_complex)
738  end do
739  end do
740  end do
Here is the caller graph for this function:

◆ create_initial_transposition_description()

type(pencil_transposition) function pencil_fft_mod::create_initial_transposition_description ( type(model_state_type), intent(inout)  current_state)
private

Creates an initial transposition representation of the Z pencil that MONC is normally decomposed in. This is then fed into the create transposition procedure which will generate transpositions to other pencils.

Parameters
current_stateThe current model state

Definition at line 614 of file pencilfft.F90.

614  type(model_state_type), intent(inout) :: current_state
615 
616  create_initial_transposition_description%dim=z_index
617  create_initial_transposition_description%process_decomposition_layout=current_state%parallel%dim_sizes
618  create_initial_transposition_description%my_process_location=current_state%parallel%my_coords
619  create_initial_transposition_description%my_pencil_size=current_state%local_grid%size
Here is the caller graph for this function:

◆ create_transposition()

type(pencil_transposition) function pencil_fft_mod::create_transposition ( type(global_grid_type), intent(inout)  global_grid,
type(pencil_transposition), intent(in)  existing_transposition,
integer, intent(in)  new_pencil_dim,
integer, dimension(:), intent(in)  process_dim_sizes,
integer, intent(in)  direction,
integer, dimension(:), intent(in)  extended_dimensions 
)
private

Creates a specific pencil transposition description. It is maybe more a decomposition description, but the main complexity comes from the transposition from existing decomposition to new decomposition so therefore it is called transposition. The new pencil decomposition depends not only on the dimension to split on, but also the existing pencil decomposition. The new decomposed dimension (i.e. the existing pencil dimension) * other local dimensions is used as the sending size, receiving though requires knowledge about the data size on the source process so others will send their pencil dimension size to this process.

Parameters
new_pencil_dimThe dimension to use as the new pencil decomposition
existing_pencil_dimThe dimension used in the current decomposition
existing_pencil_process_layoutNumber of processes per dimension for the current decomposition
existing_my_locationThe current processes block location per dimension for the current decomposition
existing_pencil_sizePencil size per dimension for the current decomposition
process_dim_sizesSizes of the pencil dimension from other processes that is used to calculate receive count
directionWhether we are transposing forwards or backwards, backwards is just an inverse
extended_dimensionsThe dimensions that this process extends from n to (n/2+1)*2 (i.e. result of fft complex->real)

Definition at line 217 of file pencilfft.F90.

217  type(global_grid_type), intent(inout) :: global_grid
218  type(pencil_transposition), intent(in) :: existing_transposition
219  integer, dimension(:), intent(in) :: process_dim_sizes
220  integer, intent(in) :: new_pencil_dim, direction, extended_dimensions(:)
221 
222  create_transposition%process_decomposition_layout=determine_pencil_process_dimensions(&
223  new_pencil_dim, existing_transposition%dim, existing_transposition%process_decomposition_layout)
224 
225  create_transposition%my_process_location=determine_my_pencil_location(new_pencil_dim, &
226  existing_transposition%dim, existing_transposition%my_process_location)
227 
228  create_transposition%my_pencil_size=determine_pencil_size(new_pencil_dim, create_transposition%process_decomposition_layout,&
229  create_transposition%my_process_location, existing_transposition, global_grid, extended_dimensions)
230 
231  allocate(create_transposition%send_dims(3, create_transposition%process_decomposition_layout(existing_transposition%dim)), &
232  create_transposition%recv_dims(3, create_transposition%process_decomposition_layout(existing_transposition%dim)))
233  if (direction == forward) then
234  call determine_my_process_sizes_per_dim(existing_transposition%dim, &
235  existing_transposition%my_pencil_size, create_transposition%process_decomposition_layout, &
236  global_grid, extended_dimensions, create_transposition%send_dims)
237  call determine_matching_process_dimensions(new_pencil_dim, existing_transposition%dim, process_dim_sizes, &
238  create_transposition%my_pencil_size, create_transposition%process_decomposition_layout, create_transposition%recv_dims)
239  else
240  call determine_my_process_sizes_per_dim(new_pencil_dim, create_transposition%my_pencil_size, &
241  existing_transposition%process_decomposition_layout, global_grid, extended_dimensions, create_transposition%recv_dims)
242  call determine_matching_process_dimensions(existing_transposition%dim, new_pencil_dim, process_dim_sizes, &
243  existing_transposition%my_pencil_size, existing_transposition%process_decomposition_layout, &
244  create_transposition%send_dims)
245  end if
246 
247  allocate(create_transposition%send_sizes(size(create_transposition%send_dims, 2)), &
248  create_transposition%send_offsets(size(create_transposition%send_sizes)), &
249  create_transposition%recv_sizes(size(create_transposition%recv_dims, 2)), &
250  create_transposition%recv_offsets(size(create_transposition%recv_sizes)))
251 
252  call concatenate_dimension_sizes(create_transposition%send_dims, create_transposition%send_sizes)
253  call determine_offsets_from_size(create_transposition%send_sizes, create_transposition%send_offsets)
254 
255  call concatenate_dimension_sizes(create_transposition%recv_dims, create_transposition%recv_sizes)
256  call determine_offsets_from_size(create_transposition%recv_sizes, create_transposition%recv_offsets)
257  create_transposition%dim=new_pencil_dim
Here is the call graph for this function:
Here is the caller graph for this function:

◆ deduce_my_global_start()

integer function pencil_fft_mod::deduce_my_global_start ( type(model_state_type), intent(inout)  current_state,
integer, intent(in)  dimension 
)
private

Determines my global start coordinate in Fourier space. This is required for cos y and cos x calculation which is fed into the tridiagonal solver. After the forward FFTs, each process has ((n/2+1)/p+r) * 2 elements, where p is the number of processes and r is the uneven process remainder (1 or 0 depending on p). Therefore some processes will have t elements, and some t-2 elements to feed into the solver.

Parameters
current_stateThe current model state
dimensionThe dimension that we are calculating this for (Y or X)
Returns
My global start in Fourier space

Definition at line 752 of file pencilfft.F90.

752  type(model_state_type), intent(inout) :: current_state
753  integer, intent(in) :: dimension
754 
755  integer complex_size, distributed_size, remainder, larger_nums, smaller_nums
756 
757  complex_size=(current_state%global_grid%size(dimension)/2+1)*2
758  distributed_size=complex_size / current_state%parallel%dim_sizes(dimension)
759  remainder=complex_size - distributed_size * current_state%parallel%dim_sizes(dimension)
760  larger_nums=min(remainder, current_state%parallel%my_coords(dimension))
761  smaller_nums=current_state%parallel%my_coords(dimension)-remainder
762  deduce_my_global_start=((distributed_size+1)*larger_nums + merge(distributed_size*smaller_nums, 0, smaller_nums .gt. 0)) + 1
Here is the caller graph for this function:

◆ determine_matching_process_dimensions()

subroutine pencil_fft_mod::determine_matching_process_dimensions ( integer, intent(in)  new_pencil_dim,
integer, intent(in)  existing_pencil_dim,
integer, dimension(:), intent(in)  proc_sizes,
integer, dimension(:), intent(in)  my_pencil_size,
integer, dimension(:), intent(in)  pencil_processes_per_dim,
integer, dimension(:,:), intent(inout)  specific_sizes_per_dim 
)
private

Determines the sizes per dimension on the matching process either to receive from (forward transposition) or send to (backwards transposition) each source process. Not only does this depend on the my pencil sizes, but it also depends on the amount of data that the source process has to send over.

Parameters
new_pencil_dimThe dimension for the new pencil decomposition
existing_pencil_dimDimension for the existing pencil decomposition
proc_sizesSize of dimension on the source processes (index in array corresponds to source PID)
my_pencil_sizeMy (new) pencil size per dimension
pencil_processes_per_dimThe process layout per dimension

Definition at line 594 of file pencilfft.F90.

594  integer, intent(in) :: new_pencil_dim, existing_pencil_dim, proc_sizes(:), my_pencil_size(:), pencil_processes_per_dim(:)
595  integer, dimension(:,:), intent(inout) :: specific_sizes_per_dim
596 
597  integer :: i, j
598 
599  do i=1,pencil_processes_per_dim(existing_pencil_dim)
600  do j=1,3
601  if (j==new_pencil_dim) then
602  specific_sizes_per_dim(j, i)=proc_sizes(i)
603  else
604  specific_sizes_per_dim(j, i)=my_pencil_size(j)
605  end if
606  end do
607  end do
Here is the caller graph for this function:

◆ determine_my_pencil_location()

integer function, dimension(3) pencil_fft_mod::determine_my_pencil_location ( integer, intent(in)  new_pencil_dim,
integer, intent(in)  existing_pencil_dim,
integer, dimension(3), intent(in)  existing_locations 
)
private

Determines my location for each dimension in the new pencil decomposition. I.e. which block I am operating on.

Parameters
new_pencil_dimNew pencil decomposition dimension
existing_pencil_dimCurrent pencil dimension
existing_locationsLocation for the current decomposition

Definition at line 554 of file pencilfft.F90.

554  integer, intent(in) :: new_pencil_dim, existing_pencil_dim, existing_locations(3)
555  integer :: determine_my_pencil_location(3)
556 
557  integer :: i
558 
559  do i=1,3
560  if (i == new_pencil_dim) then
561  determine_my_pencil_location(i)=1
562  else if (i == existing_pencil_dim) then
563  determine_my_pencil_location(i)=existing_locations(new_pencil_dim)
564  else
565  determine_my_pencil_location(i)=existing_locations(i)
566  end if
567  end do
Here is the caller graph for this function:

◆ determine_my_process_sizes_per_dim()

subroutine pencil_fft_mod::determine_my_process_sizes_per_dim ( integer, intent(in)  existing_pencil_dim,
integer, dimension(:), intent(in)  existing_pencil_size,
integer, dimension(:), intent(in)  new_pencil_procs_per_dim,
type(global_grid_type), intent(inout)  global_grid,
integer, dimension(:), intent(in)  extended_dimensions,
integer, dimension(:,:), intent(inout)  specific_sizes_per_dim 
)
private

Determines the number of elements to on my process per dimension which either need to be sent to (forwards transformation) or received from (backwards) each target process (in the row or column) This depends on the existing pencil decomposition, as effectively we are breaking that contigulity and decomposing it into n blocks in that dimension now (provided by new_pencil_procs_per_dim)

Parameters
existing_pencil_dimThe pencil dimension that we are transforming from
existing_pencil_sizeExisting pencil decomposition sizes per dimension
new_pencil_procs_per_dimFor the target decomposition the number of processes per dimension
global_gridDescription of the global grid which we use for sizing information
extended_dimensionsList of dimensions where we extend from n to n+2 (i.e. result of FFT complex-> real transformation)

Definition at line 490 of file pencilfft.F90.

490  integer, intent(in) :: existing_pencil_dim, existing_pencil_size(:), new_pencil_procs_per_dim(:), extended_dimensions(:)
491  type(global_grid_type), intent(inout) :: global_grid
492  integer, dimension(:,:), intent(inout) :: specific_sizes_per_dim
493 
494  integer :: i, split_size, split_remainder, j, s
495 
496  do i=1,3
497  if (i == existing_pencil_dim) then
498  s=global_grid%size(i)
499  if (is_extended_dimension(i, extended_dimensions)) s=s+2
500  split_size = s / new_pencil_procs_per_dim(i)
501  split_remainder = s - split_size * new_pencil_procs_per_dim(i)
502  do j=1,new_pencil_procs_per_dim(existing_pencil_dim)
503  specific_sizes_per_dim(i,j)=merge(split_size+1, split_size, j .le. split_remainder)
504  end do
505  else
506  specific_sizes_per_dim(i,:) = existing_pencil_size(i)
507  end if
508  end do
Here is the call graph for this function:
Here is the caller graph for this function:

◆ determine_offsets_from_size()

subroutine pencil_fft_mod::determine_offsets_from_size ( integer, dimension(:), intent(in)  source_sizes,
integer, dimension(:), intent(inout)  determined_offsets 
)
private

Simple helper function to deduce send or receive offsets from the sizes.

Parameters
source_sizesSizes that we are using to build the offsets

Definition at line 514 of file pencilfft.F90.

514  integer, intent(in) :: source_sizes(:)
515  integer, dimension(:), intent(inout) :: determined_offsets
516 
517  integer :: i
518 
519  determined_offsets(1)=0
520  do i=2,size(source_sizes)
521  determined_offsets(i)=determined_offsets(i-1)+source_sizes(i-1)
522  end do
Here is the caller graph for this function:

◆ determine_pencil_process_dimensions()

integer function, dimension(3) pencil_fft_mod::determine_pencil_process_dimensions ( integer, intent(in)  new_pencil_dim,
integer, intent(in)  existing_pencil_dim,
integer, dimension(3), intent(in)  existing_pencil_procs 
)
private

Determines the number of processes in each dimension for the target decomposition. This depends heavily on the existing decomposition, as we basically contiguise our pencil dimension and decompose the existing pencil dimension. The third dimension remains unchanged.

Parameters
new_pencil_dimNew pencil dimension
existing_pencil_dimCurrent decomposition pencil dimension
existing_pencil_procsCurrent decomposition process layout

Definition at line 532 of file pencilfft.F90.

532  integer, intent(in) :: new_pencil_dim, existing_pencil_dim, existing_pencil_procs(3)
533  integer :: determine_pencil_process_dimensions(3)
534 
535  integer :: i
536 
537  do i=1,3
538  if (i == new_pencil_dim) then
539  determine_pencil_process_dimensions(i)=1
540  else if (i == existing_pencil_dim) then
541  determine_pencil_process_dimensions(i)=existing_pencil_procs(new_pencil_dim)
542  else
543  determine_pencil_process_dimensions(i)=existing_pencil_procs(i)
544  end if
545  end do
Here is the caller graph for this function:

◆ determine_pencil_size()

integer function, dimension(3) pencil_fft_mod::determine_pencil_size ( integer, intent(in)  new_pencil_dim,
integer, dimension(3), intent(in)  pencil_process_layout,
integer, dimension(3), intent(in)  my_pencil_location,
type(pencil_transposition), intent(in)  existing_transposition,
type(global_grid_type), intent(inout)  global_grid,
integer, dimension(:), intent(in)  extended_dimensions 
)
private

Deduces the size of my (local) pencil based upon the new decomposition. This depends heavily on the current pencil decomposition, the new pencil dimension is the global size, the existing pencil dimension becomes decomposed based on the number of processes in that dimension. The third dimension remains unchanged.

Parameters
new_pencil_dimDimension for the new pencil decomposition
pencil_process_layoutThe processes per dimension layout for the new decomposition
my_pencil_locationMy location in the block layout
existing_pencil_dimCurrent decomposition dimension
existing_pencil_sizeCurrent decomposition sizes
global_gridDescription of the global grid which we use for sizing information
extended_dimensionsList of dimensions where we extend from n to n+2 (i.e. result of FFT complex-> real transformation)

Definition at line 634 of file pencilfft.F90.

634 
635  type(pencil_transposition), intent(in) :: existing_transposition
636  integer, intent(in) :: new_pencil_dim, pencil_process_layout(3), my_pencil_location(3), extended_dimensions(:)
637  type(global_grid_type), intent(inout) :: global_grid
638  integer :: determine_pencil_size(3)
639 
640  integer :: i, split_size, split_remainder, s
641 
642  do i=1,3
643  if (i == new_pencil_dim) then
644  if (is_extended_dimension(i, extended_dimensions)) then
645  ! If complex and Y dim then /2+1 for the global size
646  determine_pencil_size(i)=(global_grid%size(new_pencil_dim)/2+1)*2
647  else
648  determine_pencil_size(i)=global_grid%size(new_pencil_dim)
649  end if
650  else if (i == existing_transposition%dim) then
651  s=global_grid%size(i)
652  ! If complex and Y dim then use s/2+1 for the size to split
653  if (is_extended_dimension(i, extended_dimensions)) s=(s/2+1)*2
654  split_size=s/pencil_process_layout(i)
655  split_remainder=s - split_size * pencil_process_layout(i)
656  determine_pencil_size(i)=merge(split_size+1, split_size, my_pencil_location(i)+1 .le. split_remainder)
657  else
658  determine_pencil_size(i)=existing_transposition%my_pencil_size(i)
659  end if
660  end do
Here is the call graph for this function:
Here is the caller graph for this function:

◆ finalise_pencil_fft()

subroutine, public pencil_fft_mod::finalise_pencil_fft ( integer, intent(in)  monc_communicator)

Cleans up allocated buffer memory.

Definition at line 93 of file pencilfft.F90.

93  integer, intent(in) :: monc_communicator
94  integer :: ierr, i
95 
96  do i=1,size(fftw_plan_initialised)
97  if (fftw_plan_initialised(i)) then
98  call fftw_destroy_plan(fftw_plan(i))
99  end if
100  end do
101 
102  if (dim_y_comm .ne. mpi_comm_self .and. dim_y_comm .ne. monc_communicator) call mpi_comm_free(dim_y_comm, ierr)
103  if (dim_x_comm .ne. mpi_comm_self .and. dim_x_comm .ne. monc_communicator) call mpi_comm_free(dim_x_comm, ierr)
104  deallocate(buffer1, buffer2, real_buffer1, real_buffer2, real_buffer3, fft_in_y_buffer , fft_in_x_buffer)

◆ initialise_buffers()

subroutine pencil_fft_mod::initialise_buffers ( )
private

Initialises memory for the buffers used in the FFT.

Definition at line 153 of file pencilfft.F90.

153  allocate(buffer1(y_from_z_transposition%my_pencil_size(y_index)/2+1, y_from_z_transposition%my_pencil_size(x_index), &
154  y_from_z_transposition%my_pencil_size(z_index)), &
155  real_buffer1((y_from_z_transposition%my_pencil_size(y_index)/2+1)*2, y_from_z_transposition%my_pencil_size(x_index), &
156  y_from_z_transposition%my_pencil_size(z_index)), &
157  buffer2(x_from_y_transposition%my_pencil_size(x_index)/2+1, x_from_y_transposition%my_pencil_size(z_index), &
158  x_from_y_transposition%my_pencil_size(y_index)), &
159  real_buffer2((x_from_y_transposition%my_pencil_size(x_index)/2+1)*2, x_from_y_transposition%my_pencil_size(z_index), &
160  x_from_y_transposition%my_pencil_size(y_index)), &
161  fft_in_y_buffer(y_from_z_transposition%my_pencil_size(y_index), y_from_z_transposition%my_pencil_size(x_index), &
162  y_from_z_transposition%my_pencil_size(z_index)), &
163  fft_in_x_buffer(x_from_y_transposition%my_pencil_size(x_index), x_from_y_transposition%my_pencil_size(z_index), &
164  x_from_y_transposition%my_pencil_size(y_index)), &
165  real_buffer3(y_from_x_transposition%my_pencil_size(y_index), y_from_x_transposition%my_pencil_size(x_index), &
166  y_from_x_transposition%my_pencil_size(z_index)))
Here is the caller graph for this function:

◆ initialise_pencil_fft()

integer function, dimension(3), public pencil_fft_mod::initialise_pencil_fft ( type(model_state_type), intent(inout)  current_state,
integer, intent(out)  my_y_start,
integer, intent(out)  my_x_start 
)

Initialises the pencil FFT functionality, this will create the transposition structures needed.

Parameters
current_stateThe current model state
my_y_startMy global start in fourier space for Y
my_x_startMy global start in fourier space for X
Returns
Size of local dimensions in fourier space for this process

Definition at line 52 of file pencilfft.F90.

52  type(model_state_type), intent(inout) :: current_state
53  integer, intent(out) :: my_y_start, my_x_start
54  integer :: initialise_pencil_fft(3)
55 
56  integer :: ierr, y_distinct_sizes(current_state%parallel%dim_sizes(y_index)), &
57  x_distinct_sizes(current_state%parallel%dim_sizes(x_index))
58 
59  my_y_start=deduce_my_global_start(current_state, y_index)
60  my_x_start=deduce_my_global_start(current_state, x_index)
61 
62  if (current_state%parallel%dim_sizes(y_index) .gt. 1 .and. current_state%parallel%dim_sizes(x_index) .gt. 1) then
63  call mpi_cart_sub(current_state%parallel%neighbour_comm, (/1,0/), dim_y_comm, ierr)
64  call mpi_cart_sub(current_state%parallel%neighbour_comm, (/0,1/), dim_x_comm, ierr)
65 
66  call mpi_allgather(current_state%local_grid%size(y_index), 1, mpi_int, y_distinct_sizes, 1, mpi_int, dim_y_comm, ierr)
67  call mpi_allgather(current_state%local_grid%size(x_index), 1, mpi_int, x_distinct_sizes, 1, mpi_int, dim_x_comm, ierr)
68  else if (current_state%parallel%dim_sizes(y_index) .gt. 1) then
69  dim_y_comm=current_state%parallel%monc_communicator
70  dim_x_comm=mpi_comm_self
71  call mpi_allgather(current_state%local_grid%size(y_index), 1, mpi_int, y_distinct_sizes, 1, mpi_int, dim_y_comm, ierr)
72  x_distinct_sizes=current_state%local_grid%size(x_index)
73  else if (current_state%parallel%dim_sizes(x_index) .gt. 1) then
74  dim_y_comm=mpi_comm_self
75  dim_x_comm=current_state%parallel%monc_communicator
76  y_distinct_sizes=current_state%local_grid%size(y_index)
77  call mpi_allgather(current_state%local_grid%size(x_index), 1, mpi_int, x_distinct_sizes, 1, mpi_int, dim_x_comm, ierr)
78  else
79  dim_y_comm=mpi_comm_self
80  dim_x_comm=mpi_comm_self
81  y_distinct_sizes=current_state%local_grid%size(y_index)
82  x_distinct_sizes=current_state%local_grid%size(x_index)
83  end if
84 
85  call initialise_transpositions(current_state, y_distinct_sizes, x_distinct_sizes)
86  call initialise_buffers()
87 
88  initialise_pencil_fft=z_from_y_transposition%my_pencil_size
Here is the call graph for this function:

◆ initialise_transpositions()

subroutine pencil_fft_mod::initialise_transpositions ( type(model_state_type), intent(inout)  current_state,
integer, dimension(:), intent(in)  y_distinct_sizes,
integer, dimension(:), intent(in)  x_distinct_sizes 
)
private

Initialises the pencil transpositions, from a pencil in one dimension to that in another.

Parameters
current_stateThe current model state
y_distinct_sizesY sizes per process
x_distinct_sizesX sizes per process

Definition at line 174 of file pencilfft.F90.

174  type(model_state_type), intent(inout) :: current_state
175  integer, dimension(:), intent(in) :: y_distinct_sizes, x_distinct_sizes
176 
177  type(pencil_transposition) :: z_pencil
178 
179  z_pencil=create_initial_transposition_description(current_state)
180 
181  ! Transpositions
182  y_from_z_transposition=create_transposition(current_state%global_grid, z_pencil, y_index, y_distinct_sizes, &
183  forward, (/ -1 /))
184  x_from_y_transposition=create_transposition(current_state%global_grid, y_from_z_transposition, x_index, &
185  x_distinct_sizes, forward, (/ y_index /))
186  y_from_x_transposition=create_transposition(current_state%global_grid, x_from_y_transposition, y_index, &
187  normal_to_extended_process_dim_sizes(x_distinct_sizes), backward, (/ y_index, x_index /))
188  z_from_y_transposition=create_transposition(current_state%global_grid, y_from_x_transposition, z_index, &
189  normal_to_extended_process_dim_sizes(y_distinct_sizes), backward, (/ y_index, x_index /))
190 
191  y_from_z_2_transposition=create_transposition(current_state%global_grid, z_from_y_transposition, y_index, &
192  normal_to_extended_process_dim_sizes(y_distinct_sizes), forward, (/ y_index, x_index /))
193  x_from_y_2_transposition=create_transposition(current_state%global_grid, y_from_z_2_transposition, x_index, &
194  normal_to_extended_process_dim_sizes(x_distinct_sizes), forward, (/ y_index, x_index /))
195  y_from_x_2_transposition=create_transposition(current_state%global_grid, x_from_y_2_transposition, y_index, &
196  x_distinct_sizes, backward, (/ y_index /))
197  z_from_y_2_transposition=create_transposition(current_state%global_grid, y_from_x_2_transposition, z_index, &
198  y_distinct_sizes, backward, (/ -1 /))
Here is the call graph for this function:
Here is the caller graph for this function:

◆ is_extended_dimension()

logical function pencil_fft_mod::is_extended_dimension ( integer, intent(in)  dimension,
integer, dimension(:), intent(in)  extended_dimensions 
)
private

Determines whether or not the specific dimension is in the list of extended dimensions.

Parameters
dimensionThe dimension to test for
extended_dimensionsArray of dimensions that will be searched
Returns
Whether the dimension is found in the array

Definition at line 668 of file pencilfft.F90.

668  integer, intent(in) :: dimension, extended_dimensions(:)
669 
670  integer :: i
671  do i=1,size(extended_dimensions)
672  if (extended_dimensions(i) == dimension) then
673  is_extended_dimension=.true.
674  return
675  end if
676  end do
677  is_extended_dimension=.false.
Here is the caller graph for this function:

◆ normal_to_extended_process_dim_sizes()

integer function, dimension(size(process_dim_sizes)) pencil_fft_mod::normal_to_extended_process_dim_sizes ( integer, dimension(:), intent(in)  process_dim_sizes)
private

Transforms real process dimension sizes into their real after FFT complex->real transformation. The way this works is that it goes from n to (n/2+1)*2 numbers which is distributed amongst the processes deterministically.

Parameters
process_dim_sizesReal process dimension sizes
Returns
The extended process dimension sizes

Definition at line 685 of file pencilfft.F90.

685  integer, dimension(:), intent(in) :: process_dim_sizes
686  integer, dimension(size(process_dim_sizes)) :: normal_to_extended_process_dim_sizes
687 
688  integer :: temp_total, split_size, remainder
689 
690  temp_total=(sum(process_dim_sizes) /2 + 1) * 2
691  split_size=temp_total/size(process_dim_sizes)
692  remainder=temp_total - split_size*size(process_dim_sizes)
693 
694  normal_to_extended_process_dim_sizes=split_size
695  normal_to_extended_process_dim_sizes(1:remainder)=split_size+1
Here is the caller graph for this function:

◆ perform_backwards_3dfft()

subroutine, public pencil_fft_mod::perform_backwards_3dfft ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:,:,:), intent(in)  source_data,
real(kind=default_precision), dimension(:,:,:), intent(out)  target_data 
)

Performs a backwards 3D FFT and currently results in target data which is the X, Z, Y oriented pencil Note that the source_data here takes no account for the halo, it is up to caller to exclude this. This does no FFT in Z, but transposes to Y, does FFT in Y, then transposes to X and performs an FFT in that dimension. Pencil decomposition is used which has already been set up.

Parameters
current_stateThe current model state
source_dataThe source real data to in the frequency domain
target_dataTime domain complex representation of the frequency domain source

Definition at line 138 of file pencilfft.F90.

138  type(model_state_type), target, intent(inout) :: current_state
139  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(in) :: source_data
140  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(out) :: target_data
141 
142  call transpose_to_pencil(y_from_z_2_transposition, (/z_index, y_index, x_index/), dim_y_comm, forward, &
143  source_data, real_buffer3)
144  call transpose_to_pencil(x_from_y_2_transposition, (/y_index, x_index, z_index/), dim_x_comm, forward, &
145  real_buffer3, real_buffer2)
146 
147  call transpose_and_backward_fft_in_x(current_state, real_buffer2, buffer2, real_buffer1)
148  call transpose_and_backward_fft_in_y(current_state, real_buffer1, buffer1, target_data)
Here is the call graph for this function:

◆ perform_c2r_fft()

subroutine pencil_fft_mod::perform_c2r_fft ( complex(c_double_complex), dimension(:,:,:), intent(inout), pointer, contiguous  source_data,
real(kind=default_precision), dimension(:,:,:), intent(inout), pointer, contiguous  transformed_data,
integer, intent(in)  row_size,
integer, intent(in)  num_rows,
integer, intent(in)  plan_id 
)
private

Performs the complex to real (backwards) FFT.

Parameters
source_dataSource (complex) data in the frequency domain
transformed_dataResulting real data in the time domain
row_sizeNumber of elements for each FFT
num_rowsThe number of FFTs to perform on the next data elements in the source_data
plan_idId number of the plan that tracks whether we need to create it or can reuse the existing one

Definition at line 450 of file pencilfft.F90.

450  complex(C_DOUBLE_COMPLEX), dimension(:,:,:), contiguous, pointer, intent(inout) :: source_data
451  real(kind=DEFAULT_PRECISION), dimension(:,:,:), contiguous, pointer, intent(inout) :: transformed_data
452  integer, intent(in) :: row_size, num_rows, plan_id
453 
454  if (.not. fftw_plan_initialised(plan_id)) then
455  ! n is the size of the FFT (in real, not complex->real coords.) There are row_size/2+1 between entries for the input
456  ! (complex) data and row_size between entries for the output data
457  fftw_plan(plan_id) = fftw_plan_many_dft_c2r(1, (/row_size/), num_rows, source_data, (/row_size/2+1/), 1, row_size/2+1, &
458  transformed_data, (/row_size/), 1, row_size, fftw_estimate)
459  fftw_plan_initialised(plan_id)=.true.
460  end if
461  call fftw_execute_dft_c2r(fftw_plan(plan_id), source_data, transformed_data)
Here is the caller graph for this function:

◆ perform_forward_3dfft()

subroutine, public pencil_fft_mod::perform_forward_3dfft ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:,:,:), intent(inout)  source_data,
real(kind=default_precision), dimension(:,:,:), intent(out)  target_data 
)

Performs a forward 3D FFT and currently results in target data which is the X, Z, Y oriented pencil Note that the source_data here takes no account for the halo, it is up to caller to exclude this. This does no FFT in Z, but transposes to Y, does FFT in Y, then transposes to X and performs an FFT in that dimension. Pencil decomposition is used which has already been set up.

Parameters
current_stateThe current model state
source_dataThe source real data to in the time domain
target_dataFrequency domain real representation of the time domain source which is allocated here

Definition at line 115 of file pencilfft.F90.

115  type(model_state_type), target, intent(inout) :: current_state
116  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(inout) :: source_data
117  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(out) :: target_data
118 
119  call transpose_and_forward_fft_in_y(current_state, source_data, buffer1, real_buffer1)
120  real_buffer1=real_buffer1/current_state%global_grid%size(y_index)
121  call transpose_and_forward_fft_in_x(current_state, real_buffer1, buffer2, real_buffer2)
122  real_buffer2=real_buffer2/current_state%global_grid%size(x_index)
123 
124  call transpose_to_pencil(y_from_x_transposition, (/x_index, z_index, y_index/), dim_x_comm, backward, &
125  real_buffer2, real_buffer3)
126  call transpose_to_pencil(z_from_y_transposition, (/y_index, x_index, z_index/), dim_y_comm, backward, &
127  real_buffer3, target_data)
Here is the call graph for this function:

◆ perform_r2c_fft()

subroutine pencil_fft_mod::perform_r2c_fft ( real(kind=default_precision), dimension(:,:,:), intent(inout), pointer, contiguous  source_data,
complex(c_double_complex), dimension(:,:,:), intent(inout), pointer, contiguous  transformed_data,
integer, intent(in)  row_size,
integer, intent(in)  num_rows,
integer, intent(in)  plan_id 
)
private

Actually performs a forward real to complex FFT.

Parameters
source_dataSource (real) data in the time domain
transformed_dataResulting complex data in the frequency domain
row_sizeNumber of elements for each FFT
num_rowsThe number of FFTs to perform on the next data elements in the source_data
plan_idId number of the plan that tracks whether we need to create it or can reuse the existing one

Definition at line 431 of file pencilfft.F90.

431  real(kind=DEFAULT_PRECISION), dimension(:,:,:), contiguous, pointer, intent(inout) :: source_data
432  complex(C_DOUBLE_COMPLEX), dimension(:,:,:), contiguous, pointer, intent(inout) :: transformed_data
433  integer, intent(in) :: row_size, num_rows, plan_id
434 
435  if (.not. fftw_plan_initialised(plan_id)) then
436  fftw_plan(plan_id) = fftw_plan_many_dft_r2c(1, (/row_size/), num_rows, source_data, (/row_size/), 1, row_size, &
437  transformed_data, (/row_size/), 1, row_size/2+1, fftw_estimate)
438  fftw_plan_initialised(plan_id)=.true.
439  end if
440  call fftw_execute_dft_r2c(fftw_plan(plan_id), source_data, transformed_data)
Here is the caller graph for this function:

◆ rearrange_data_for_sending()

subroutine pencil_fft_mod::rearrange_data_for_sending ( real(kind=default_precision), dimension(:,:,:), intent(in)  real_source,
real(kind=default_precision), dimension(:,:,:), intent(out)  real_target 
)
private

Rearranges data for sending, transposing a,b,c into c,b,a . This is done as alltoall splits on dimension c so to go from one pencil to another we assume here that a is the existing pencil as it is contiguous.

Parameters
real_sourceSource data to transpose from
real_targetTarget data to transpose to

Definition at line 469 of file pencilfft.F90.

469  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(in) :: real_source
470  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(out) :: real_target
471 
472  integer :: i
473 
474  do i=1, size(real_source,2)
475  real_target(:,i,:)=transpose(real_source(:,i,:))
476  end do
Here is the caller graph for this function:

◆ transpose_and_backward_fft_in_x()

subroutine pencil_fft_mod::transpose_and_backward_fft_in_x ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:,:,:), intent(inout)  source_data,
complex(c_double_complex), dimension(:,:,:), intent(out), pointer, contiguous  buffer,
real(kind=default_precision), dimension(:,:,:), intent(out)  real_buffer 
)
private

Performs the backwards FFT in X and then transposes to Y pencil. The FFT requires complex numbers which are converted to real, so the this real to complex operation is performed first. If n is the logical size of the FFT row, then the input size is n+2, complex number size is n/2+1 and we get n reals out.

Parameters
current_stateThe current model state
source_dataInput buffer, X pencil oriented x,z,y
bufferComplex buffer which is fed into the FFT
real_bufferOutput buffer, Y pencil, oriented y,x,z

Definition at line 289 of file pencilfft.F90.

289  type(model_state_type), target, intent(inout) :: current_state
290  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(inout) :: source_data
291  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(out) :: real_buffer
292  complex(C_DOUBLE_COMPLEX), dimension(:,:,:), contiguous, pointer, intent(out) :: buffer
293 
294  call convert_real_to_complex(source_data, buffer)
295  call perform_c2r_fft(buffer, fft_in_x_buffer, x_from_y_2_transposition%my_pencil_size(x_index)-2, &
296  x_from_y_2_transposition%my_pencil_size(y_index) * x_from_y_2_transposition%my_pencil_size(z_index), 2)
297 
298  ! Transpose globally from X pencil to Y pencil
299  call transpose_to_pencil(y_from_x_2_transposition, (/x_index, z_index, y_index/), dim_x_comm, backward, &
300  fft_in_x_buffer, real_buffer)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ transpose_and_backward_fft_in_y()

subroutine pencil_fft_mod::transpose_and_backward_fft_in_y ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:,:,:), intent(inout)  source_data,
complex(c_double_complex), dimension(:,:,:), intent(out), pointer, contiguous  buffer,
real(kind=default_precision), dimension(:,:,:), intent(out)  real_buffer 
)
private

Performs the backwards FFT in Y and then transposes to Z pencil. The FFT requires complex numbers which are converted to real, so the this real to complex operation is performed first. If n is the logical size of the FFT row, then the input size is n+2, complex number size is n/2+1 and we get n reals out.

Parameters
current_stateThe current model state
source_dataInput buffer, Y pencil oriented y,x,z
bufferComplex buffer which is fed into the FFT
real_bufferOutput buffer, Z pencil, oriented z,y,x

Definition at line 332 of file pencilfft.F90.

332  type(model_state_type), target, intent(inout) :: current_state
333  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(inout) :: source_data
334  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(out) :: real_buffer
335  complex(C_DOUBLE_COMPLEX), dimension(:,:,:), contiguous, pointer, intent(out) :: buffer
336 
337  call convert_real_to_complex(source_data, buffer)
338 
339  call perform_c2r_fft(buffer, fft_in_y_buffer, y_from_x_2_transposition%my_pencil_size(y_index)-2, &
340  y_from_x_2_transposition%my_pencil_size(x_index) * y_from_x_2_transposition%my_pencil_size(z_index), 4)
341 
342  ! Go from global Y pencil to global Z pencil
343  call transpose_to_pencil(z_from_y_2_transposition, (/y_index, x_index, z_index/), dim_y_comm, backward, &
344  fft_in_y_buffer, real_buffer)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ transpose_and_forward_fft_in_x()

subroutine pencil_fft_mod::transpose_and_forward_fft_in_x ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:,:,:), intent(inout)  buffer1,
complex(c_double_complex), dimension(:,:,:), intent(out), pointer, contiguous  buffer2,
real(kind=default_precision), dimension(:,:,:), intent(inout)  buffer3 
)
private

Performs the transposition and forward FFT in the x dimension. After the FFT the complex space is converted back into real numbers. The X size is (n/2+1)*2 due to this transformation.

Parameters
current_stateThe current model state
buffer1Input buffer, Y pencil after the Y dimension FFT oriented y,x,z
bufferComplex buffer which results from the FFT
buffer2Output buffer, X pencil after this X FFT, oriented x,z,y

Definition at line 310 of file pencilfft.F90.

310  type(model_state_type), target, intent(inout) :: current_state
311  complex(C_DOUBLE_COMPLEX), dimension(:,:,:), contiguous, pointer, intent(out) :: buffer2
312  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(inout) :: buffer1, buffer3
313 
314  ! Go from global Y pencil to global X pencil
315  call transpose_to_pencil(x_from_y_transposition, (/y_index, x_index, z_index/), dim_x_comm, forward, &
316  buffer1, fft_in_x_buffer)
317 
318  call perform_r2c_fft(fft_in_x_buffer, buffer2, x_from_y_transposition%my_pencil_size(x_index), &
319  x_from_y_transposition%my_pencil_size(y_index) * x_from_y_transposition%my_pencil_size(z_index), 3)
320 
321  call convert_complex_to_real(buffer2, buffer3)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ transpose_and_forward_fft_in_y()

subroutine pencil_fft_mod::transpose_and_forward_fft_in_y ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:,:,:), intent(inout)  source_data,
complex(c_double_complex), dimension(:,:,:), intent(out), pointer, contiguous  buffer,
real(kind=default_precision), dimension(:,:,:), intent(out)  real_buffer 
)
private

Performs the transposition and forward FFT in the y dimension then converts back to real numbers. The Y size is (n/2+1)*2 due to the complex to real transformation after the FFT.

Parameters
current_stateThe current model state
source_dataInput buffer, Z pencil oriented z,y,x
bufferComplex buffer which the FFT writes into
real_bufferOutput buffer, Y pencil, oriented y,x,z

Definition at line 267 of file pencilfft.F90.

267  type(model_state_type), target, intent(inout) :: current_state
268  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(inout) :: source_data
269  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(out) :: real_buffer
270  complex(C_DOUBLE_COMPLEX), dimension(:,:,:), contiguous, pointer, intent(out) :: buffer
271 
272  ! Transpose globally from Z pencil to Y pencil
273  call transpose_to_pencil(y_from_z_transposition, (/z_index, y_index, x_index/), dim_y_comm, forward, &
274  source_data, fft_in_y_buffer)
275 
276  call perform_r2c_fft(fft_in_y_buffer, buffer, y_from_z_transposition%my_pencil_size(y_index), &
277  y_from_z_transposition%my_pencil_size(x_index) * y_from_z_transposition%my_pencil_size(z_index), 1)
278  call convert_complex_to_real(buffer, real_buffer)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ transpose_to_pencil()

subroutine pencil_fft_mod::transpose_to_pencil ( type(pencil_transposition), intent(in)  transposition_description,
integer, dimension(3), intent(in)  source_dims,
integer, intent(in)  communicator,
integer, intent(in)  direction,
real(kind=default_precision), dimension(:,:,:), intent(in)  source_data,
real(kind=default_precision), dimension(:,:,:), intent(out)  target_data 
)
private

Transposes globally to a new pencil decomposition. This goes from the source dimensions a,b,c to b,c,a (forwards) or c,a,b (backwards.) It requires multiple steps, first the local data is transposed to c,b,a regardless of direction. then it is communicated via alltoall, each process then assembles its own b,c,a or c,a,b data via contiguising across blocks as the data layout is nonlinear.

Parameters
transposition_descriptionDescription of the transposition
source_dimsDimensions of the current pencil that we wish to transpose from, will go from abc to bca
communicatorThe MPI communicator associated with the group of processes who will swap data
directionWhether this is going forwards or backwards, it makes a difference to the data arrangement
source_dataSource data (abc)
target_dataTarget data (bca)

Definition at line 358 of file pencilfft.F90.

358  type(pencil_transposition), intent(in) :: transposition_description
359  integer, intent(in) :: source_dims(3), communicator, direction
360  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(in) :: source_data
361  real(kind=DEFAULT_PRECISION), dimension(:,:,:), intent(out) :: target_data
362 
363  integer :: ierr
364  real(kind=DEFAULT_PRECISION), dimension(:,:,:), allocatable :: real_temp
365  real(kind=DEFAULT_PRECISION), dimension(:), allocatable :: real_temp2
366 
367 
368  allocate(real_temp(size(source_data,3), size(source_data,2), size(source_data,1)), &
369  real_temp2(product(transposition_description%my_pencil_size)+1))
370 
371  call rearrange_data_for_sending(real_source=source_data, real_target=real_temp)
372 
373  call mpi_alltoallv(real_temp, transposition_description%send_sizes, transposition_description%send_offsets, &
374  precision_type, real_temp2, transposition_description%recv_sizes, transposition_description%recv_offsets, &
375  precision_type, communicator, ierr)
376  call contiguise_data(transposition_description, (/source_dims(3), source_dims(2), source_dims(1)/), direction, &
377  source_real_buffer=real_temp2, target_real_buffer=target_data)
378  deallocate(real_temp, real_temp2)
Here is the call graph for this function:
Here is the caller graph for this function:

Variable Documentation

◆ backward

integer, parameter pencil_fft_mod::backward =2
private

Transposition directions.

Definition at line 28 of file pencilfft.F90.

◆ buffer1

complex(c_double_complex), dimension(:,:,:), pointer, contiguous pencil_fft_mod::buffer1
private

Definition at line 37 of file pencilfft.F90.

37  complex(C_DOUBLE_COMPLEX), dimension(:,:,:), contiguous, pointer :: buffer1, buffer2

◆ buffer2

complex(c_double_complex), dimension(:,:,:), pointer, contiguous pencil_fft_mod::buffer2
private

Definition at line 37 of file pencilfft.F90.

◆ dim_x_comm

integer pencil_fft_mod::dim_x_comm
private

Communicators for each dimension.

Definition at line 29 of file pencilfft.F90.

◆ dim_y_comm

integer pencil_fft_mod::dim_y_comm
private

Definition at line 29 of file pencilfft.F90.

29  integer :: dim_y_comm, dim_x_comm

◆ fft_in_x_buffer

real(kind=default_precision), dimension(:,:,:), pointer, contiguous pencil_fft_mod::fft_in_x_buffer
private

Definition at line 35 of file pencilfft.F90.

◆ fft_in_y_buffer

real(kind=default_precision), dimension(:,:,:), pointer, contiguous pencil_fft_mod::fft_in_y_buffer
private

Definition at line 35 of file pencilfft.F90.

◆ fftw_plan

type(c_ptr), dimension(4) pencil_fft_mod::fftw_plan
private

Definition at line 40 of file pencilfft.F90.

40  type(c_ptr) :: fftw_plan(4)

◆ fftw_plan_initialised

logical, dimension(4) pencil_fft_mod::fftw_plan_initialised =.false.
private

Definition at line 41 of file pencilfft.F90.

41  logical :: fftw_plan_initialised(4)=.false.

◆ forward

integer, parameter pencil_fft_mod::forward =1
private

Definition at line 28 of file pencilfft.F90.

28  integer, parameter :: forward=1, backward=2

◆ real_buffer1

real(kind=default_precision), dimension(:,:,:), pointer, contiguous pencil_fft_mod::real_buffer1
private

Definition at line 35 of file pencilfft.F90.

35  real(kind=DEFAULT_PRECISION), dimension(:,:,:), contiguous, pointer :: real_buffer1, real_buffer2, real_buffer3, &
36  fft_in_y_buffer , fft_in_x_buffer

◆ real_buffer2

real(kind=default_precision), dimension(:,:,:), pointer, contiguous pencil_fft_mod::real_buffer2
private

Definition at line 35 of file pencilfft.F90.

◆ real_buffer3

real(kind=default_precision), dimension(:,:,:), pointer, contiguous pencil_fft_mod::real_buffer3
private

Definition at line 35 of file pencilfft.F90.

◆ x_from_y_2_transposition

type(pencil_transposition) pencil_fft_mod::x_from_y_2_transposition
private

Definition at line 31 of file pencilfft.F90.

◆ x_from_y_transposition

type(pencil_transposition) pencil_fft_mod::x_from_y_transposition
private

Definition at line 31 of file pencilfft.F90.

◆ y_from_x_2_transposition

type(pencil_transposition) pencil_fft_mod::y_from_x_2_transposition
private

Definition at line 31 of file pencilfft.F90.

◆ y_from_x_transposition

type(pencil_transposition) pencil_fft_mod::y_from_x_transposition
private

Definition at line 31 of file pencilfft.F90.

◆ y_from_z_2_transposition

type(pencil_transposition) pencil_fft_mod::y_from_z_2_transposition
private

Definition at line 31 of file pencilfft.F90.

◆ y_from_z_transposition

type(pencil_transposition) pencil_fft_mod::y_from_z_transposition
private

Definition at line 31 of file pencilfft.F90.

31  type(pencil_transposition) :: y_from_z_transposition, x_from_y_transposition, y_from_x_transposition, z_from_y_transposition, &
32  y_from_z_2_transposition, x_from_y_2_transposition, y_from_x_2_transposition, z_from_y_2_transposition

◆ z_from_y_2_transposition

type(pencil_transposition) pencil_fft_mod::z_from_y_2_transposition
private

Definition at line 31 of file pencilfft.F90.

◆ z_from_y_transposition

type(pencil_transposition) pencil_fft_mod::z_from_y_transposition
private

Definition at line 31 of file pencilfft.F90.