draw_markers: send the markers to their processes also without sorting boxes

src/struphy/pic/base.py

@@ -1146,12 +1146,14 @@ def draw_markers(
 
         All needed parameters can be set in the parameter file, see :ref:`params_yml`.
 
-        An initial sorting will be performed if sort is given as True (default) and sorting_params were given to the init.
+        On several processes, the markers are always sent to the process of their domain after drawing
+        (:meth:`mpi_sort_markers`). They are also sorted into boxes if sort is True (default) and sorting
+        boxes were set up.
 
         Parameters
         ----------
         sort : Bool
-            Wether to sort the particules in boxes after initial drawing (only if sorting params were passed)
+            Whether to sort the markers into boxes after initial drawing (only if sorting boxes were set up).
         """
 
         import scipy.special as sp
@@ -1405,10 +1407,15 @@ def draw_markers(
             assert xp.all(~self.holes[:n_mks_load_loc])
             assert xp.all(self.holes[n_mks_load_loc:])
 
-        if self._initialized_sorting and sort:
-            logger.info("\nSorting the markers after initial draw")
-            if self.mpi_comm is not None:
-                self.mpi_sort_markers()
+        # send the markers to the process of their domain: needed on several processes whether or not there are
+        # sorting boxes, since accumulations before the first push (e.g. an initial Poisson solve) assume it (#764)
+        sort_boxes = self._initialized_sorting and sort
+        if self.mpi_comm is not None and (self.mpi_size > 1 or sort_boxes):
+            logger.info("\nSending the markers to their processes after initial draw")
+            self.mpi_sort_markers()
+
+        if sort_boxes:
+            logger.info("\nSorting the markers into boxes after initial draw")
             self.do_sort()
             logger.info("Done.")
 

src/struphy/pic/tests/test_draw_parallel.py

@@ -177,6 +177,56 @@ def test_marker_ids(loading, Np=1000):
     assert n_at_origin == 0
 
 
+@pytest.mark.parametrize("loading", ["pseudo_random", "sobol_standard", "sobol_antithetic"])
+@pytest.mark.parametrize("sort", [True, False])
+def test_draw_without_boxes_on_rank_domain(loading, sort, Np=1000):
+    """Without sorting boxes, every marker is on the domain of its process right after drawing (#764),
+    for ``sort=True`` and ``sort=False`` (which only concerns the boxes), and no marker is lost."""
+
+    import cunumpy as xp
+    from maybempi import MPI
+
+    from struphy import LoadingParameters, domains
+    from struphy.feec.psydac_derham import Derham
+    from struphy.io.options import DerhamOptions
+    from struphy.pic.particles import Particles6D
+    from struphy.topology.grids import TensorProductGrid
+
+    comm = MPI.COMM_WORLD
+
+    grid = TensorProductGrid(num_elements=[8, 5, 3])
+    derham = Derham(grid, DerhamOptions(degree=[1, 1, 1]), comm=comm)
+    domain_decomp = (derham.domain_array, derham.domain_decomposition.nprocs)
+
+    loading_params = LoadingParameters(
+        Np=Np,
+        seed=1234,
+        loading=loading,
+        moments=(0.0, 0.0, 0.0, 1.0, 1.0, 1.0),
+    )
+
+    # no sorting_params: no sorting boxes
+    particles = Particles6D(
+        comm_world=comm,
+        domain_decomp=domain_decomp,
+        loading_params=loading_params,
+        domain=domains.Cuboid(),
+    )
+    assert particles.sorting_boxes is None
+
+    particles.draw_markers(sort=sort)
+
+    positions = particles.markers_wo_holes[:, :3]
+    on_rank = xp.all(particles._is_on_domain_of_rank(positions, particles.mpi_rank), axis=1)
+    n_off_rank = int(xp.sum(~on_rank))
+    n_markers = positions.shape[0]
+    if particles.mpi_size > 1:
+        n_markers = particles.mpi_comm.allreduce(n_markers)
+
+    assert n_off_rank == 0, f"rank {particles.mpi_rank}: {n_off_rank} markers off the process domain after drawing"
+    assert n_markers == Np
+
+
 if __name__ == "__main__":
     # test_draw([8, 9, 10], [2, 3, 4], [False, False, True], ['Cuboid', {
     #     'l1': 1., 'r1': 2., 'l2': 10., 'r2': 20., 'l3': 100., 'r3': 200.}])