Decompose across panels (#220)

Co-authored-by: Hacka Fett <christophermaynard@users.noreply.github.com>

Author: thomasmelvin

CONTRIBUTORS.md

@@ -20,3 +20,4 @@
 | t00sa            | Sam Clarke-Green       | Met Office  | 2026-01-27 |
 | MetBenjaminWent  | Benjamin Went          | Met Office  | 2026-01-30 |
 | jcsmeto          | James Cunningham-Smith | Met Office  | 2026-02-06 |
+| thomasmelvin     | Thomas Melvin          | Met Office  | 2026-01-15 |

components/driver/source/mesh/runtime_partition_mod.f90

@@ -76,15 +76,15 @@ subroutine get_partition_strategy( mesh_selection, total_ranks, partitioner_ptr
       call log_event( "Using serial cubed sphere partitioner", &
                       log_level_debug )
 
-    else if (mod(total_ranks,6) == 0) then
+    else if (mod(total_ranks,3) == 0 .or. mod(total_ranks,2) == 0) then
       ! Paralled run job
       partitioner_ptr => partitioner_cubedsphere
       call log_event( "Using parallel cubed sphere partitioner", &
                       log_level_debug )
 
     else
-      call log_event( "Total number of processors must be 1 (serial) "// &
-                      "or a multiple of 6 for a cubed-sphere domain.",   &
+      call log_event( "Total number of processors must be 1 (serial) "//    &
+                      "or a multiple of 2 or 3 for a cubed-sphere domain.", &
                       log_level_error )
     end if
 

documentation/source/how_to_use_it/technical_articles/lfric_distmem_impl.rst

@@ -177,8 +177,12 @@ supplied to the partitioner through the use of a function passed in via
 a function pointer.
 
 For efficiency reasons, the algorithm for partitioning a cubed-sphere
-mesh is a specific algorithm that will only work for a cubed-sphere mesh
-(each panel of the cubed-sphere is split into rectangular partitions).
+mesh is a specific algorithm that will only work for a cubed-sphere mesh.
+By default each panel of the cubed-sphere is split into rectangular partitions,
+alternatively if the 'custom' decomposition type is used then sets of 2 or 3 
+cubed-sphere panels can be grouped together before these sets are decomposed into
+a number of rectangular panels. These restrictions mean that the total number of
+partitions for cubed-sphere meshes needs to be a multiple of 2 or 3.
 This specific partitioner could be replaced with a more general
 partitioner for use on fully unstructured meshes.
 

infrastructure/source/mesh/global_mesh_mod.F90

@@ -1019,79 +1019,104 @@ end function get_equatorial_latitude
   !>
   !> @note For a cubed -sphere mesh, this will only return correct cell IDs if
   !>       the offset cell remains on same the cubed-sphere "face" as the
-  !>       start cell.
+  !>       start cell unless check_orientation is specified.
   !>
   !> @param[in] cell_number ID of the anchor element.
   !> @param[in] x_cells Offset in the E/W direction.
   !> @param[in] y_cells Offset in the N/S direction.
+  !> @param[in] check_orientation Switch to check for orientation changes
+  !!                              such as when crossing panel boundaries
   !>
   !> @return cell_id ID of the cell at the given offset to the start cell.
   !>
   function get_cell_id( self, cell_number, &
-                        x_cells, y_cells ) result ( cell_id )
+                        x_cells, y_cells,  &
+                        check_orientation ) result ( cell_id )
 
   use reference_element_mod, only : W, S, E, N
 
   implicit none
 
   class(global_mesh_type), intent(in) :: self
 
-  integer(i_def), intent(in) :: cell_number
-  integer(i_def), intent(in) :: x_cells, y_cells
+  integer(i_def),           intent(in) :: cell_number
+  integer(i_def),           intent(in) :: x_cells, y_cells
+  logical(l_def), optional, intent(in) :: check_orientation
 
-  integer(i_def) :: cell_id
+  integer(i_def) :: cell_id, old_cell_id
 
-  integer(i_def) :: index_x, dist_x
-  integer(i_def) :: index_y, dist_y
-  integer(i_def) :: i
+  integer(i_def) :: x_index, y_index, x_dist, y_dist, i, j
+  integer(i_def) :: opposite(4), rotate(4)
 
-  cell_id = cell_number
+  logical(l_def) :: check
 
-  ! Determine march along local x-axis
-  if (x_cells > 0) then
-    index_x = E
-    dist_x  = x_cells
-  else if (x_cells < 0) then
-    index_x = W
-    dist_x  = abs(x_cells)
+  opposite = (/ E, N, W, S /)
+  rotate = (/ S, E, N, W /)
+  if ( present( check_orientation ) ) then
+    check = check_orientation
   else
-    index_x = W
-    dist_x  = 0
+    check = .false.
   end if
 
-  ! Determine march along local y-axis
-  if (y_cells > 0) then
-    index_y = N
-    dist_y  = y_cells
-  else if (y_cells < 0) then
-    index_y = S
-    dist_y  = abs(y_cells)
+  cell_id=cell_number
+
+  if (x_cells >= 0 )then
+    x_index = E
+    x_dist = x_cells
+  else
+    x_index = W
+    x_dist = abs(x_cells)
+  end if
+  if (y_cells >= 0 )then
+    y_index = N
+    y_dist = y_cells
   else
-    index_y = S
-    dist_y  = 0
+    y_index = S
+    y_dist = abs(y_cells)
   end if
 
-  !========================================
-  ! March from anchor along local x/y axes.
-  !========================================
-  do i=1, dist_x
+  ! x_dist cells in the x-direction
+  do i = 1, x_dist
+    old_cell_id = cell_id
+    cell_id = self%cell_next_2d(x_index,cell_id)
     if (cell_id == self%void_cell) then
       ! The current cell is not on the domain
       write(log_scratch_space,'(A)') &
           'No adjacent cell, the current cell is not on mesh domain'
       call log_event(log_scratch_space, LOG_LEVEL_ERROR)
     end if
-    cell_id = self%cell_next_2d(index_x, cell_id)
+    ! Check if we've changed direction
+    if ( check .and. self%cell_next_2d(opposite(x_index),cell_id) /= old_cell_id ) then
+      ! We have changed direction so we need to find the correct
+      ! index and reset
+      do j = 1,4
+        if ( self%cell_next_2d(opposite(j),cell_id) == old_cell_id ) x_index = j
+      end do
+    end if
   end do
-
-  do i=1, dist_y
+  ! y_dist cells in the y-direction
+  ! Since the direction may have changed we need to recompute
+  y_index = rotate(x_index)
+  if ( y_cells < 0 ) y_index = opposite(y_index)
+
+  ! y_index and y_dist
+  do i = 1,y_dist
+    old_cell_id = cell_id
+    cell_id = self%cell_next_2d(y_index,cell_id)
     if (cell_id == self%void_cell) then
       ! The current cell is not on the domain
       write(log_scratch_space,'(A)') &
           'No adjacent cell, the current cell is not on mesh domain'
       call log_event(log_scratch_space, LOG_LEVEL_ERROR)
     end if
-    cell_id = self%cell_next_2d(index_y, cell_id)
+    ! Check if we've changed direction
+    if ( check .and.  self%cell_next_2d(opposite(y_index),cell_id) /= old_cell_id ) then
+      ! We have changed direction so we need to find the correct
+      ! index and reset
+      do j = 1,4
+        if ( self%cell_next_2d(opposite(j),cell_id) == old_cell_id ) y_index = j
+      end do
+    end if
   end do
 
   end function get_cell_id