submodule(mesh_cart_m) mesh_cart_communicate_s !! Submodule related with MPI communication of mesh use mpi_mapping_auxiliaries_m, only : getdata_fwdbwdplane implicit none contains module subroutine exchange_mesh_mpi(self, comm, step, mesh_received) class(mesh_cart_t), intent(in) :: self integer, intent(in) :: comm integer, intent(in) :: step type(mesh_cart_t), intent(out) :: mesh_received integer :: ns, nr, szn, yperiodic_int, extend_beyond_wall_int integer, allocatable, dimension(:) :: iwrk_send, iwrk_recv real(FP), allocatable, dimension(:) :: rwrk_send, rwrk_recv ! Communicate integer and logical metadata if (self%yperiodic) then yperiodic_int = 1 else yperiodic_int = 0 endif if (self%extend_beyond_wall) then extend_beyond_wall_int = 1 else extend_beyond_wall_int = 0 endif ns = 18 allocate(iwrk_send(ns)) iwrk_send = [self%lvl, & self%lvst, & self%n_points, & self%nx_f, & self%ny_f, & self%min_i, & self%max_i, & self%min_j, & self%max_j, & self%size_neighbor, & self%size_ghost_layer, & self%n_points_inner, & self%n_points_boundary, & self%n_points_ghost, & self%n_points_red, & self%n_points_black, & yperiodic_int, & extend_beyond_wall_int] call getdata_fwdbwdplane(comm, step, ns, iwrk_send, nr, iwrk_recv) mesh_received%lvl = iwrk_recv(1) mesh_received%lvst = iwrk_recv(2) mesh_received%n_points = iwrk_recv(3) mesh_received%nx_f = iwrk_recv(4) mesh_received%ny_f = iwrk_recv(5) mesh_received%min_i = iwrk_recv(6) mesh_received%max_i = iwrk_recv(7) mesh_received%min_j = iwrk_recv(8) mesh_received%max_j = iwrk_recv(9) mesh_received%size_neighbor = iwrk_recv(10) mesh_received%size_ghost_layer = iwrk_recv(11) mesh_received%n_points_inner = iwrk_recv(12) mesh_received%n_points_boundary = iwrk_recv(13) mesh_received%n_points_ghost = iwrk_recv(14) mesh_received%n_points_red = iwrk_recv(15) mesh_received%n_points_black = iwrk_recv(16) if (iwrk_recv(17) == 0) then mesh_received%yperiodic = .false. else mesh_received%yperiodic = .true. endif if (iwrk_recv(18) == 0) then mesh_received%extend_beyond_wall = .false. else mesh_received%extend_beyond_wall = .true. endif deallocate(iwrk_recv) deallocate(iwrk_send) ! Communicate real metadata ns = 5 allocate(rwrk_send(ns)) rwrk_send = [self%spacing_f, & self%spacing_c, & self%phi, & self%xmin, & self%ymin] call getdata_fwdbwdplane(comm, step, ns, rwrk_send, nr, rwrk_recv) mesh_received%spacing_f = rwrk_recv(1) mesh_received%spacing_c = rwrk_recv(2) mesh_received%phi = rwrk_recv(3) mesh_received%xmin = rwrk_recv(4) mesh_received%ymin = rwrk_recv(5) deallocate(rwrk_recv) deallocate(rwrk_send) ! Communicate 1D datasets ns = self%n_points call getdata_fwdbwdplane(comm, step, ns, self%cart_i, & nr, mesh_received%cart_i) call getdata_fwdbwdplane(comm, step, ns, self%cart_j, & nr, mesh_received%cart_j) call getdata_fwdbwdplane(comm, step, ns, self%pinfo, & nr, mesh_received%pinfo) call getdata_fwdbwdplane(comm, step, ns, self%district, & nr, mesh_received%district) call getdata_fwdbwdplane(comm, step, ns, self%redblack_indices, & nr, mesh_received%redblack_indices) ns = self%n_points_inner call getdata_fwdbwdplane(comm, step, ns, self%inner_indices, & nr, mesh_received%inner_indices) ns = self%n_points_boundary call getdata_fwdbwdplane(comm, step, ns, self%boundary_indices, & nr, mesh_received%boundary_indices) ns = self%n_points_ghost call getdata_fwdbwdplane(comm, step, ns, self%ghost_indices, & nr, mesh_received%ghost_indices) ! Communicate index_neighbor (multi-dimensional) ns = self%n_points * (2 * self%size_neighbor + 1)**2 allocate(iwrk_send(ns), source=0) ! iwrk_send is initialised here with source=0 values explixitly ! as it is accessed directly below in shape function. If it was only ! allocated, problems could arise due to compiler optimisations iwrk_send = reshape(self%index_neighbor, shape(iwrk_send)) call getdata_fwdbwdplane(comm, step, ns, iwrk_send, & nr, iwrk_recv) szn = mesh_received%size_neighbor allocate(mesh_received%index_neighbor(-szn:szn, -szn:szn, & mesh_received%n_points), & source=0) mesh_received%index_neighbor = & reshape(iwrk_recv, shape(mesh_received%index_neighbor)) deallocate(iwrk_recv) deallocate(iwrk_send) ! For cart_shaped_l array check firstly if allocated on ! mesh to be fetched and then build explicitly ! (A direct communication might get lock, if allocated ! on some ranks, but not on others) ns = 1 allocate(iwrk_send(ns)) if (allocated(self%cart_shaped_l)) then iwrk_send(1) = 1 else iwrk_send(1) = 0 endif call getdata_fwdbwdplane(comm, step, ns, iwrk_send, & nr, iwrk_recv) if (iwrk_recv(1) == 1) then call mesh_received%build_shaped_cart_arrays() endif deallocate(iwrk_recv) deallocate(iwrk_send) end subroutine end submodule