From abddc5fd716ed91ae2a747c839331f9e31d7cdb6 Mon Sep 17 00:00:00 2001 From: nmannall Date: Tue, 14 Jan 2025 15:25:01 +0000 Subject: [PATCH] Update singleuse commands for parallel model build --- gprMax/contexts.py | 5 +- gprMax/grid/fdtd_grid.py | 33 +++++--- gprMax/grid/mpi_grid.py | 108 +++++++++++++++++++++----- gprMax/model.py | 5 +- gprMax/mpi_model.py | 22 ++---- gprMax/scene.py | 24 +++--- gprMax/user_inputs.py | 18 ++--- gprMax/user_objects/cmds_multiuse.py | 3 +- gprMax/user_objects/cmds_singleuse.py | 43 +++++----- 9 files changed, 167 insertions(+), 94 deletions(-) diff --git a/gprMax/contexts.py b/gprMax/contexts.py index 9f9dc421..11668d56 100644 --- a/gprMax/contexts.py +++ b/gprMax/contexts.py @@ -201,9 +201,8 @@ class MPIContext(Context): if not model_config.reuse_geometry(): model = self._create_model() - if model.is_coordinator(): - scene = self._get_scene(model_num) - scene.create_internal_objects(model) + scene = self._get_scene(model_num) + scene.create_internal_objects(model) model.build() diff --git a/gprMax/grid/fdtd_grid.py b/gprMax/grid/fdtd_grid.py index 0813da0f..331c3ebf 100644 --- a/gprMax/grid/fdtd_grid.py +++ b/gprMax/grid/fdtd_grid.py @@ -21,7 +21,7 @@ import itertools import logging import sys from collections import OrderedDict -from typing import Any, Iterable, List, Tuple, Union +from typing import Any, Iterable, List, Literal, Tuple, Union import numpy as np import numpy.typing as npt @@ -98,7 +98,8 @@ class FDTDGrid: # itself does not matter, however, if must be the same from model to # model otherwise the numerical precision from adding the PML # corrections will be different. - self.pmls["thickness"] = OrderedDict((key, 10) for key in PML.boundaryIDs) + self.pmls["thickness"] = OrderedDict() + self.set_pml_thickness(10) # Materials used by this grid self.materials: List[Material] = [] @@ -140,6 +141,18 @@ class FDTDGrid: def dz(self, value: float): self.dl[2] = value + def set_pml_thickness(self, thickness: Union[int, Tuple[int, int, int, int, int, int]]): + if isinstance(thickness, int) or len(thickness) == 1: + for key in PML.boundaryIDs: + self.pmls["thickness"][key] = int(thickness) + elif len(thickness) == 6: + self.pmls["thickness"]["x0"] = int(thickness[0]) + self.pmls["thickness"]["y0"] = int(thickness[1]) + self.pmls["thickness"]["z0"] = int(thickness[2]) + self.pmls["thickness"]["xmax"] = int(thickness[3]) + self.pmls["thickness"]["ymax"] = int(thickness[4]) + self.pmls["thickness"]["zmax"] = int(thickness[5]) + def build(self) -> None: """Build the grid.""" @@ -435,15 +448,15 @@ class FDTDGrid: logger.exception("Receiver(s) will be stepped to a position outside the domain.") raise ValueError from e - IntPoint = Tuple[int, int, int] - FloatPoint = Tuple[float, float, float] - - def within_bounds(self, p: IntPoint): + def within_bounds(self, p: Tuple[int, int, int]) -> bool: """Check a point is within the grid. Args: p: Point to check. + Returns: + within_bounds: True if the point is within the grid bounds. + Raises: ValueError: Raised if the point is outside the grid. """ @@ -454,7 +467,9 @@ class FDTDGrid: if p[2] < 0 or p[2] > self.nz: raise ValueError("z") - def discretise_point(self, p: FloatPoint) -> IntPoint: + return True + + def discretise_point(self, p: Tuple[float, float, float]) -> Tuple[int, int, int]: """Calculate the nearest grid cell to the given point. Args: @@ -468,7 +483,7 @@ class FDTDGrid: z = round_value(float(p[2]) / self.dz) return (x, y, z) - def round_to_grid(self, p: FloatPoint) -> FloatPoint: + def round_to_grid(self, p: Tuple[float, float, float]) -> Tuple[float, float, float]: """Round the provided point to the nearest grid cell. Args: @@ -481,7 +496,7 @@ class FDTDGrid: p_r = (p[0] * self.dx, p[1] * self.dy, p[2] * self.dz) return p_r - def within_pml(self, p: IntPoint) -> bool: + def within_pml(self, p: Tuple[int, int, int]) -> bool: """Check if the provided point is within a PML. Args: diff --git a/gprMax/grid/mpi_grid.py b/gprMax/grid/mpi_grid.py index 1fc5cc26..88797a22 100644 --- a/gprMax/grid/mpi_grid.py +++ b/gprMax/grid/mpi_grid.py @@ -59,8 +59,6 @@ class MPIGrid(FDTDGrid): def __init__(self, comm: MPI.Cartcomm): self.size = np.zeros(3, dtype=np.intc) - super().__init__() - self.comm = comm self.x_comm = comm.Sub([False, True, True]) self.y_comm = comm.Sub([True, False, True]) @@ -82,6 +80,8 @@ class MPIGrid(FDTDGrid): self.send_halo_map = np.empty((3, 2), dtype=MPI.Datatype) self.recv_halo_map = np.empty((3, 2), dtype=MPI.Datatype) + super().__init__() + @property def rank(self) -> int: return self.comm.Get_rank() @@ -114,6 +114,47 @@ class MPIGrid(FDTDGrid): def nz(self, value: int): self.size[Dim.Z] = value + @property + def gx(self) -> int: + return self.global_size[Dim.X] + + @gx.setter + def gx(self, value: int): + self.global_size[Dim.X] = value + + @property + def gy(self) -> int: + return self.global_size[Dim.Y] + + @gy.setter + def gy(self, value: int): + self.global_size[Dim.Y] = value + + @property + def gz(self) -> int: + return self.global_size[Dim.Z] + + @gz.setter + def gz(self, value: int): + self.global_size[Dim.Z] = value + + def set_pml_thickness(self, thickness: Union[int, Tuple[int, int, int, int, int, int]]): + super().set_pml_thickness(thickness) + + # Set PML thickness to zero if not at the edge of the domain + if self.has_neighbour(Dim.X, Dir.NEG): + self.pmls["thickness"]["x0"] = 0 + if self.has_neighbour(Dim.X, Dir.POS): + self.pmls["thickness"]["xmax"] = 0 + if self.has_neighbour(Dim.Y, Dir.NEG): + self.pmls["thickness"]["y0"] = 0 + if self.has_neighbour(Dim.Y, Dir.POS): + self.pmls["thickness"]["ymax"] = 0 + if self.has_neighbour(Dim.Z, Dir.NEG): + self.pmls["thickness"]["z0"] = 0 + if self.has_neighbour(Dim.Z, Dir.POS): + self.pmls["thickness"]["zmax"] = 0 + def is_coordinator(self) -> bool: """Test if the current rank is the coordinator. @@ -478,20 +519,6 @@ class MPIGrid(FDTDGrid): self.scatter_snapshots() - self.pmls = self.comm.bcast(self.pmls, root=self.COORDINATOR_RANK) - if self.coords[Dim.X] != 0: - self.pmls["thickness"]["x0"] = 0 - if self.coords[Dim.X] != self.mpi_tasks[Dim.X] - 1: - self.pmls["thickness"]["xmax"] = 0 - if self.coords[Dim.Y] != 0: - self.pmls["thickness"]["y0"] = 0 - if self.coords[Dim.Y] != self.mpi_tasks[Dim.Y] - 1: - self.pmls["thickness"]["ymax"] = 0 - if self.coords[Dim.Z] != 0: - self.pmls["thickness"]["z0"] = 0 - if self.coords[Dim.Z] != self.mpi_tasks[Dim.Z] - 1: - self.pmls["thickness"]["zmax"] = 0 - if not self.is_coordinator(): # TODO: When scatter arrays properly, should initialise these to the local grid size self.initialise_geometry_arrays() @@ -666,7 +693,6 @@ class MPIGrid(FDTDGrid): ) raise ValueError - self.calculate_local_extents() self.set_halo_map() self.distribute_grid() @@ -751,6 +777,50 @@ class MPIGrid(FDTDGrid): self.upper_extent = self.lower_extent + self.size logger.debug( - f"Local grid size: {self.size}, Lower extent: {self.lower_extent}, Upper extent:" - f" {self.upper_extent}" + f"Global grid size: {self.global_size}, Local grid size: {self.size}, Lower extent:" + f" {self.lower_extent}, Upper extent: {self.upper_extent}" + ) + + def within_bounds(self, p: Tuple[int, int, int]) -> bool: + """Check a point is within the grid. + + Args: + p: Point to check. + + Returns: + within_bounds: True if the point is within the local grid + (i.e. this rank's grid). False otherwise. + + Raises: + ValueError: Raised if the point is outside the global grid. + """ + if p[0] < 0 or p[0] > self.gx: + raise ValueError("x") + if p[1] < 0 or p[1] > self.gy: + raise ValueError("y") + if p[2] < 0 or p[2] > self.gz: + raise ValueError("z") + + local_point = self.global_to_local_coordinate(np.array(p, dtype=np.intc)) + + return all(local_point >= self.negative_halo_offset) and all(local_point <= self.size) + + def within_pml(self, p: Tuple[int, int, int]) -> bool: + """Check if the provided point is within a PML. + + Args: + p: Point to check. + + Returns: + within_pml: True if the point is within a PML. + """ + local_point = self.global_to_local_coordinate(np.array(p)) + p = (local_point[0], local_point[1], local_point[2]) + + # within_pml check will only be valid if the point is also + # within the local grid + return ( + super().within_pml(p) + and all(local_point >= self.negative_halo_offset) + and all(local_point <= self.size) ) diff --git a/gprMax/model.py b/gprMax/model.py index 8d7012f5..b3481e1f 100644 --- a/gprMax/model.py +++ b/gprMax/model.py @@ -19,7 +19,7 @@ import datetime import logging import sys -from typing import List, Sequence +from typing import List, Sequence, Tuple import humanize import numpy as np @@ -170,6 +170,9 @@ class Model: return grid + def set_size(self, size: Tuple[int, int, int]): + self.nx, self.ny, self.nz = size + def build(self): """Builds the Yee cells for a model.""" diff --git a/gprMax/mpi_model.py b/gprMax/mpi_model.py index 6d073b68..72aee350 100644 --- a/gprMax/mpi_model.py +++ b/gprMax/mpi_model.py @@ -1,5 +1,5 @@ import logging -from typing import Optional +from typing import Optional, Tuple import numpy as np from mpi4py import MPI @@ -53,6 +53,11 @@ class MPIModel(Model): def is_coordinator(self): return self.rank == 0 + def set_size(self, size: Tuple[int, int, int]): + super().set_size(size) + + self.G.calculate_local_extents() + def build_geometry(self): self._broadcast_model() @@ -61,23 +66,8 @@ class MPIModel(Model): self._filter_geometry_objects() def _broadcast_model(self): - self.title = self.comm.bcast(self.title) - - self.nx = self.comm.bcast(self.nx) - self.ny = self.comm.bcast(self.ny) - self.nz = self.comm.bcast(self.nz) - - self.comm.Bcast(self.dl) - self.dt = self.comm.bcast(self.dt) - - self.iterations = self.comm.bcast(self.iterations) - - self.srcsteps = self.comm.bcast(self.srcsteps) - self.rxsteps = self.comm.bcast(self.rxsteps) - model_config = config.get_model_config() model_config.materials["maxpoles"] = self.comm.bcast(model_config.materials["maxpoles"]) - model_config.ompthreads = self.comm.bcast(model_config.ompthreads) def _filter_geometry_objects(self): objects = self.comm.bcast(self.geometryobjects) diff --git a/gprMax/scene.py b/gprMax/scene.py index 5f592812..6c9ff952 100644 --- a/gprMax/scene.py +++ b/gprMax/scene.py @@ -21,6 +21,7 @@ from typing import List, Sequence from gprMax.grid.fdtd_grid import FDTDGrid from gprMax.materials import create_built_in_materials from gprMax.model import Model +from gprMax.mpi_model import MPIModel from gprMax.subgrids.user_objects import SubGridBase as SubGridUserBase from gprMax.user_objects.cmds_geometry.add_grass import AddGrass from gprMax.user_objects.cmds_geometry.add_surface_roughness import AddSurfaceRoughness @@ -174,15 +175,20 @@ class Scene: # Process commands that can only have a single instance self.process_single_use_objects(model) - # Process multiple commands - self.process_multi_use_objects(model) + if ( + isinstance(model, MPIModel) + and model.is_coordinator() + or not isinstance(model, MPIModel) + ): + # Process multiple commands + self.process_multi_use_objects(model) - # Initialise geometry arrays for main and subgrids - for grid in [model.G] + model.subgrids: - grid.initialise_geometry_arrays() + # Initialise geometry arrays for main and subgrids + for grid in [model.G] + model.subgrids: + grid.initialise_geometry_arrays() - # Process the main grid geometry commands - self.process_geometry_objects(self.geometry_objects, model.G) + # Process the main grid geometry commands + self.process_geometry_objects(self.geometry_objects, model.G) - # Process all the commands for subgrids - self.process_subgrid_objects(model) + # Process all the commands for subgrids + self.process_subgrid_objects(model) diff --git a/gprMax/user_inputs.py b/gprMax/user_inputs.py index b72b2bc8..566d50a8 100644 --- a/gprMax/user_inputs.py +++ b/gprMax/user_inputs.py @@ -18,7 +18,7 @@ from __future__ import annotations import logging -from typing import Generic +from typing import Generic, Tuple import numpy as np from typing_extensions import TypeVar @@ -48,10 +48,13 @@ class UserInput(Generic[GridType]): def __init__(self, grid: GridType): self.grid = grid - def point_within_bounds(self, p, cmd_str, name): + def point_within_bounds(self, p, cmd_str, name, ignore_error=False) -> bool: try: - self.grid.within_bounds(p) + return self.grid.within_bounds(p) except ValueError as err: + if ignore_error: + return False + v = ["x", "y", "z"] # Discretisation dl = getattr(self.grid, f"d{err.args[0]}") @@ -64,17 +67,10 @@ class UserInput(Generic[GridType]): logger.exception(s) raise - def check_point_within_bounds(self, p) -> bool: - try: - self.grid.within_bounds(p) - return True - except ValueError: - return False - def grid_upper_bound(self) -> list[int]: return [self.grid.nx, self.grid.ny, self.grid.nz] - def discretise_point(self, p): + def discretise_point(self, p: Tuple[float, float, float]) -> Tuple[int, int, int]: """Gets the index of a continuous point with the grid.""" rv = np.vectorize(round_value) return rv(p / self.grid.dl) diff --git a/gprMax/user_objects/cmds_multiuse.py b/gprMax/user_objects/cmds_multiuse.py index fc34a347..6eedc95a 100644 --- a/gprMax/user_objects/cmds_multiuse.py +++ b/gprMax/user_objects/cmds_multiuse.py @@ -18,7 +18,6 @@ import inspect import logging -from abc import ABC, abstractmethod from os import PathLike from pathlib import Path from typing import Optional, Union @@ -1165,7 +1164,7 @@ class Snapshot(GridUserObject): while any(p2 < upper_bound): # Ideally extend p2 up to the correct upper bound. This will # not change the snapshot output. - if uip.check_point_within_bounds(upper_bound): + if uip.point_within_bounds(upper_bound, "", "", ignore_error=True): p2 = upper_bound p2_continuous = uip.descretised_to_continuous(p2) logger.warning( diff --git a/gprMax/user_objects/cmds_singleuse.py b/gprMax/user_objects/cmds_singleuse.py index 8ac2a6ba..5932b8e3 100644 --- a/gprMax/user_objects/cmds_singleuse.py +++ b/gprMax/user_objects/cmds_singleuse.py @@ -122,12 +122,11 @@ class Domain(ModelUserObject): def build(self, model: Model): uip = self._create_uip(model.G) - model.nx, model.ny, model.nz = uip.discretise_point(self.domain_size) - # TODO: Remove when distribute full build for MPI - if isinstance(model.G, MPIGrid): - model.G.nx = model.nx - model.G.ny = model.ny - model.G.nz = model.nz + + discretised_domain_size = uip.discretise_point(self.domain_size) + + # TODO: Fix type hinting + model.set_size(discretised_domain_size) if model.nx == 0 or model.ny == 0 or model.nz == 0: raise ValueError(f"{self} requires at least one cell in every dimension") @@ -327,7 +326,7 @@ class PMLFormulation(ModelUserObject): def build(self, model: Model): if self.formulation not in PML.formulations: - logger.exception(f"{self} requires the value to be one of {' '.join(PML.formulations)}") + raise ValueError(f"{self} requires the value to be one of {' '.join(PML.formulations)}") model.G.pmls["formulation"] = self.formulation @@ -367,27 +366,23 @@ class PMLThickness(ModelUserObject): def build(self, model: Model): grid = model.G - if isinstance(self.thickness, int) or len(self.thickness) == 1: - for key in grid.pmls["thickness"].keys(): - grid.pmls["thickness"][key] = int(self.thickness) - elif len(self.thickness) == 6: - grid.pmls["thickness"]["x0"] = int(self.thickness[0]) - grid.pmls["thickness"]["y0"] = int(self.thickness[1]) - grid.pmls["thickness"]["z0"] = int(self.thickness[2]) - grid.pmls["thickness"]["xmax"] = int(self.thickness[3]) - grid.pmls["thickness"]["ymax"] = int(self.thickness[4]) - grid.pmls["thickness"]["zmax"] = int(self.thickness[5]) - else: + if not ( + isinstance(self.thickness, int) or len(self.thickness) == 1 or len(self.thickness) == 6 + ): raise ValueError(f"{self} requires either one or six parameter(s)") + model.G.set_pml_thickness(self.thickness) + # Check each PML does not take up more than half the grid + # TODO: MPI ranks not containing a PML will not throw an error + # here. if ( - 2 * grid.pmls["thickness"]["x0"] >= grid.nx - or 2 * grid.pmls["thickness"]["y0"] >= grid.ny - or 2 * grid.pmls["thickness"]["z0"] >= grid.nz - or 2 * grid.pmls["thickness"]["xmax"] >= grid.nx - or 2 * grid.pmls["thickness"]["ymax"] >= grid.ny - or 2 * grid.pmls["thickness"]["zmax"] >= grid.nz + 2 * grid.pmls["thickness"]["x0"] >= model.nx + or 2 * grid.pmls["thickness"]["y0"] >= model.ny + or 2 * grid.pmls["thickness"]["z0"] >= model.nz + or 2 * grid.pmls["thickness"]["xmax"] >= model.nx + or 2 * grid.pmls["thickness"]["ymax"] >= model.ny + or 2 * grid.pmls["thickness"]["zmax"] >= model.nz ): raise ValueError(f"{self} has too many cells for the domain size")