Module `phi.geom`

Differentiable geometry package.

Classes:

Geometry (base type)
Box
Sphere

See the phi.geom module documentation at https://tum-pbs.github.io/PhiFlow/Geometry.html

Functions

Represent existing geometry as a signed distance function.

Args

geo: Geometry to represent as a signed distance function. Must implement Geometry.approximate_signed_distance().
bounds: Bounds of the SDF. If None will be determined from bounds of geo and rel_margin/abs_margin.
rel_margin: Relative size to pad the domain on all sides around the bounds of geo. For example, 0.1 will pad 10% of geo's size in each axis on both sides.
abs_margin: World-space size to pad the domain on all sides around the bounds of geo.
separate: Dimensions along which to unstack geo and return individual SDFs. Once created, SDFs cannot be unstacked.

Returns:

def assert_same_rank(rank1, rank2, error_message)

Tests that two objects have the same spatial rank. Objects can be of types: int, None (no check), Geometry, Shape, Tensor

def bounding_box(geometry: phi.geom._geom.Geometry | phiml.math._tensors.Tensor, reduce=<function non_batch>) ‑> phi.geom._box.Box

Builds a bounding box around geometry or a collection of points.

Args

geometry: Geometry object or Tensor of points.
reduce: Which objects to includes in each bounding box. Non-reduced dims will be part of the returned box.

Returns

Bounding Box containing only batch dims and vector.

def build_mesh(bounds: phi.geom._box.Box = None, resolution=(), obstacles: phi.geom._geom.Geometry | Dict[str, phi.geom._geom.Geometry] = None, method='quad', cell_dim: phiml.math._shape.Shape = (cellsⁱ=None), face_format: str = 'csc', max_squish: float | None = 0.5, **resolution_: int | phiml.math._tensors.Tensor | tuple | list | Any) ‑> phi.geom._mesh.Mesh

Build a mesh for a given domain, respecting obstacles.

Args

bounds: Bounds for uniform cells.
resolution: Base resolution
obstacles: Single Geometry or dict mapping boundary name to corresponding Geometry.
method: Meshing algorithm. Only quad is currently supported.
cell_dim: Dimension along which to list the cells. This should be an instance dimension.
face_format: Sparse storage format for cell connectivity.
max_squish: Smallest allowed cell size compared to the smallest regular cell.
**resolution_: For uniform grid, pass resolution as int and specify bounds. Or pass a sequence of floats for each dimension, specifying the vertex positions along each axis. This allows for variable cell stretching.

Returns

Mesh

Clips the length of a vector to the interval [min_len, max_len] while keeping the direction. Zero-vectors remain zero-vectors.

Args

vec: Tensor
min_len: Lower clipping threshold.
max_len: Upper clipping threshold.
vec_dim: Dimensions to compute the length over. By default, all channel dimensions are used to compute the vector length.
eps: Minimum vector length. Use to avoid inf gradients for zero-length vectors.

Returns

Tensor with same shape as vec.

def concat(values: Sequence[~PhiTreeNodeType], dim: phiml.math._shape.Shape | str, expand_values=False, **kwargs) ‑> ~PhiTreeNodeType

Concatenates a sequence of phiml.math.magic.Shapable objects, e.g. Tensor, along one dimension. All values must have the same spatial, instance and channel dimensions and their sizes must be equal, except for dim. Batch dimensions will be added as needed.

Args

values: Tuple or list of phiml.math.magic.Shapable, such as phiml.math.Tensor
dim: Concatenation dimension, must be present in all values. The size along dim is determined from values and can be set to undefined (None). Alternatively, a str of the form 't->name:t' can be specified, where t is on of b d i s c denoting the dimension type. This first packs all dimensions of the input into a new dim with given name and type, then concatenates the values along this dim.
expand_values: If True, will first add missing dimensions to all values, not just batch dimensions. This allows tensors with different dimensions to be concatenated. The resulting tensor will have all dimensions that are present in values.
**kwargs: Additional keyword arguments required by specific implementations. Adding spatial dimensions to fields requires the bounds: Box argument specifying the physical extent of the new dimensions. Adding batch dimensions must always work without keyword arguments.

Returns

Concatenated Tensor

Examples

>>> concat([math.zeros(batch(b=10)), math.ones(batch(b=10))], 'b')
(bᵇ=20) 0.500 ± 0.500 (0e+00...1e+00)

>>> concat([vec(x=1, y=0), vec(z=2.)], 'vector')
(x=1.000, y=0.000, z=2.000) float64

def cross(vec1: phiml.math._tensors.Tensor, vec2: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Computes the cross product of two vectors in 2D.

Args

vec1: Tensor with a single channel dimension called 'vector'
vec2: Tensor with a single channel dimension called 'vector'

Returns

Tensor

def cylinder(center: phiml.math._tensors.Tensor | float = None, radius: phiml.math._tensors.Tensor | float = None, depth: phiml.math._tensors.Tensor | float = None, rotation: phiml.math._tensors.Tensor | None = None, axis: int | str | phiml.math._tensors.Tensor = -1, variables=('center', 'radius', 'depth', 'rotation'), **center_: phiml.math._tensors.Tensor | float)

Args

center: Cylinder center as Tensor with vector dimension. The spatial dimension order should be specified in the vector dimension via item names. Can be left empty to specify dimensions via kwargs.
radius: Cylinder radius as float or Tensor.
depth: Cylinder length as float or Tensor.
rotation: Rotation angle(s) or rotation matrix.
axis: The cylinder is aligned along this axis, perturbed by rotation. Specified either as the dim along which the cylinder is aligned or as a vector.
variables: Which properties of the cylinder are variable, i.e. traced and optimizable. All by default.
**center_: Specifies center when the center argument is not given. Center position by dimension, e.g. x=0.5, y=0.2.

def embed(geometry: phi.geom._geom.Geometry, projected_dims: phiml.math._shape.Shape | str | tuple | list | None) ‑> phi.geom._geom.Geometry

Adds fake spatial dimensions to a geometry. The geometry value will be constant along the added dimensions, as if it had infinite length in these directions.

Dimensions that are already present with geometry are ignored.

Args

geometry: Geometry
projected_dims: Additional dimensions

Returns

Geometry with spatial rank geometry.spatial_rank + projected_dims.rank.

def enclosing_grid(*geometries: phi.geom._geom.Geometry, voxel_count: int, rel_margin=0.0, abs_margin=0.0) ‑> phi.geom._grid.UniformGrid

Constructs a UniformGrid which fully encloses the geometries. The grid voxels are chosen to have approximately the same size along each axis.

Args

*geometries: Geometry objects which should lie within the grid.
voxel_count: Approximate number of total voxels.
rel_margin: Relative margin, i.e. empty space on each side as a fraction of the bounding box size of geometries.
abs_margin: Absolute margin, i.e. empty space on each side.

Returns

UniformGrid

def graph(nodes: phi.geom._geom.Geometry | phiml.math._tensors.Tensor, edges: phiml.math._tensors.Tensor, boundary: Dict[str, Dict[str, slice]] = None, build_distances=True, build_bounding_distance=False) ‑> phi.geom._graph.Graph

Construct a Graph.

Args

nodes: Location Tensor or Geometry objects representing the nodes.
edges: Connectivity and edge value Tensor.
boundary: Named boundary sets.
build_distances: Whether to compute all edge lengths. This enables the properties Graph.deltas, Graph.unit_deltas, Graph.distances.
build_bounding_distance: Whether to compute the maximum edge length. This enables the property Graph.bounding_distance.

Returns

Graph

def infinite_cylinder(center=None, radius=None, inf_dim: str | phiml.math._shape.Shape | tuple | list = None, **center_) ‑> phi.geom._geom.Geometry

Creates an infinite cylinder. This is equal to embedding an n-dimensional Sphere in n+1 dimensions.

Args

center: Center coordinates without inf_dim. Alternatively use keyword arguments.
radius: Cylinder radius.
inf_dim: Dimension along which the cylinder is infinite. Use Geometry.rotated() if the direction does not align with an axis.
**center_: Alternatively specify center coordinates without inf_dim as keyword arguments.

Returns

Geometry

def intersection(*geometries: phi.geom._geom.Geometry, dim=(intersectionⁱ=None)) ‑> phi.geom._geom.Geometry

Intersection of the given geometries. A point lies inside the union if it lies within all of the geometries.

Args

*geometries: arbitrary geometries with same spatial dims. Arbitrary batch dims are allowed.
dim: Intersection dimension. This must be an instance dimension.

Returns

intersection Geometry

def invert(geometry: phi.geom._geom.Geometry)

Swaps inside and outside.

Args

geometry: Geometry to swap

Returns

New Geometry object with same surface but swapped normals

def length(obj: phi.geom._geom.Geometry | phiml.math._tensors.Tensor, epsilon=1e-05) ‑> phiml.math._tensors.Tensor

Returns the length of a vector Tensor or geometric object with a length-like property.

Args

obj:
epsilon: Minimum valid vector length. Use to avoid inf gradients for zero-length vectors. Lengths shorter than eps are set to 0.

Returns

Length as Tensor

def line_trace(geo: phi.geom._geom.Geometry, origin: phiml.math._tensors.Tensor, direction: phiml.math._tensors.Tensor, side='both', tolerance=None, max_iter=64, step_size=0.9, max_line_length=None) ‑> Tuple[phiml.math._tensors.Tensor, phiml.math._tensors.Tensor, phiml.math._tensors.Tensor, phiml.math._tensors.Tensor, phiml.math._tensors.Tensor | None]

Trace a line until it hits the surface of geo. The surface can be hit either from the outside or the inside.

Args

geo: Geometry that implements approximate_closest_surface.
origin: Line start location.
direction: Unit vector pointing in the line direction.
side: 'outside' or 'inside' or 'both'.
tolerance: Surface distance tolerance.
max_iter: Maximum number of steps per line.
step_size: Step size factor. This can be set to 1 if the signed distance values of geo are exact. For inexact SDFs, smaller step sizes prevent skipping over surfaces.

Returns

hit: Whether a surface intersection was found for the line.
distance: Distance between the line and the surface.
position: Hit location or point until which the line was traced.
normal: Surface normal at hit location
hit_index: Geometry face index at hit location

def load_gmsh(file: str, boundary_names: Sequence[str] = None, periodic: str = None, cell_dim=(cellsⁱ=None), face_format: str = 'csc')

Load an unstructured mesh from a .msh file.

This requires the package meshio to be installed.

Args

file: Path to .su2 file.
boundary_names: Boundary identifiers corresponding to the blocks in the file. If not specified, boundaries will be numbered.
periodic:
cell_dim: Dimension along which to list the cells. This should be an instance dimension.
face_format: Sparse storage format for cell connectivity.

Returns

Mesh

def load_stl(file: str, face_dim=(facesⁱ=None)) ‑> phi.geom._mesh.Mesh

Load a triangle Mesh from an STL file.

Args

file: File path to .stl file.
face_dim: Instance dim along which to list the triangles.

Returns

Mesh with spatial_rank=3 and element_rank=2.

def load_su2(file_or_mesh: str, cell_dim=(cellsⁱ=None), face_format: str = 'csc') ‑> phi.geom._mesh.Mesh

Load an unstructured mesh from a .su2 file.

This requires the package ezmesh to be installed.

Args

file_or_mesh: Path to .su2 file or ezmesh Mesh instance.
cell_dim: Dimension along which to list the cells. This should be an instance dimension.
face_format: Sparse storage format for cell connectivity.

Returns

Mesh

def mesh(vertices: phi.geom._geom.Geometry | phiml.math._tensors.Tensor, elements: phiml.math._tensors.Tensor, boundaries: str | Dict[str, List[Sequence]] | None = None, element_rank: int = None, periodic: str = None, face_format: str = 'csc', max_cell_walk: int = None)

Create a mesh from vertex positions and vertex lists.

Args

vertices: Tensor with one instance and one channel dimension vector.
elements: Lists of vertex indices as 2D tensor. The elements must be listed along an instance dimension, and the vertex indices belonging to the same polygon must be listed along a spatial dimension.
boundaries: Pass a str to assign one name to all boundary faces. For multiple boundaries, pass a dict mapping group names str to lists of faces, defined by their vertices. The last entry can be None to group all boundary faces not explicitly listed before. The boundaries dict maps boundary names to a list of edges (point pairs) in 2D and faces (3 or more points) in 3D (not yet supported).
face_format: Storage format for cell connectivity, must be one of csc, coo, csr, dense.

Returns

Mesh

def mesh_from_numpy(points: Sequence[Sequence], polygons: Sequence[Sequence], boundaries: str | Dict[str, List[Sequence]] | None = None, element_rank: int = None, periodic: str = None, cell_dim: phiml.math._shape.Shape = (cellsⁱ=None), face_format: str = 'csc') ‑> phi.geom._mesh.Mesh

Construct an unstructured mesh from vertices.

Args

points: 2D numpy array of shape (num_points, point_coord). The last dimension must have length 2 for 2D meshes and 3 for 3D meshes.
polygons: List of elements. Each polygon is defined as a sequence of point indices mapping into `points'. E.g. [(0, 1, 2)] denotes a single triangle connecting points 0, 1, and 2.
boundaries: An unstructured mesh can have multiple boundaries, each defined by a name str and a list of faces, defined by their vertices. The boundaries dict maps boundary names to a list of edges (point pairs) in 2D and faces (3 or more points) in 3D (not yet supported).
cell_dim: Dimension along which to list the cells. This should be an instance dimension.
face_format: Storage format for cell connectivity, must be one of csc, coo, csr, dense.

Returns

Mesh

def normal_from_slope(slope: phiml.math._tensors.Tensor, space: str | phiml.math._shape.Shape | Sequence[str])

Computes the normal vector of a line, plane, or hyperplane.

Args

slope: Line Slope (2D), plane slope (3D) or hyperplane slope (4+D). Must have one channel dimension listing the vector components. The vector must list all but one dimensions of space.
space: Ordered spatial dimensions as comma-separated string, sequence of names or Shape

Returns

Normal vector with the channel dimension of slope listing all dimensions of space in that order.

def normalize(obj: phiml.math._tensors.Tensor, epsilon=1e-05, allow_infinite=False, allow_zero=True)

Normalize a vector Tensor along the 'vector' dim.

Args

obj: Tensor with 'vector' dim.
epsilon: (Optional) Zero-length threshold. Vectors shorter than this length yield the unit vector (1, 0, 0, …). If not specified, the zero-vector yields NaN as it cannot be normalized.
allow_infinite: Allow infinite components in vectors. These vectors will then only points towards the infinite components.
allow_zero: Whether to return zero vectors for inputs smaller epsilon instead of a unit vector.

Returns

Tensor of the same shape as obj.

def numpy_sdf(sdf: Callable, bounds: phi.geom._box.BaseBox, center: phiml.math._tensors.Tensor = None) ‑> phi.geom._sdf.SDF

Define a SDF (signed distance function) from a NumPy function.

Args

sdf: Function mapping a location numpy.ndarray of shape (points, vector) to the corresponding SDF value (points,).
bounds: Bounds inside which the function is defined.
center: Optional center position of the object encoded via sdf.

Returns

SDF

Compresses multiple dimensions into a single dimension by concatenating the elements. Elements along the new dimensions are laid out according to the order of dims. If the order of dims differs from the current dimension order, the tensor is transposed accordingly. This function replaces the traditional reshape for these cases.

The type of the new dimension will be equal to the types of dims. If dims have varying types, the new dimension will be a batch dimension.

If none of dims exist on value, packed_dim will be added only if it is given with a definite size and value is not a primitive type.

Args

value: phiml.math.magic.Shapable, such as phiml.math.Tensor.
dims: Dimensions to be compressed in the specified order.
packed_dim: Single-dimension Shape.
pos: Index of new dimension. None for automatic, -1 for last, 0 for first.
**kwargs: Additional keyword arguments required by specific implementations. Adding spatial dimensions to fields requires the bounds: Box argument specifying the physical extent of the new dimensions. Adding batch dimensions must always work without keyword arguments.

Returns

Same type as value.

Examples

>>> pack_dims(math.zeros(spatial(x=4, y=3)), spatial, instance('points'))
(pointsⁱ=12) const 0.0

def rotate(obj: ~GeometricType, rot: float | phiml.math._tensors.Tensor | None, invert=False, pivot: str | phiml.math._tensors.Tensor = 'bounds') ‑> ~GeometricType

Rotate a vector or Geometry about the pivot.

Args

obj: n-dimensional vector Tensor or Geometry.
rot: Euler angle(s) or rotation matrix. None is interpreted as no rotation.
invert: Whether to apply the inverse rotation.
pivot: Either a point (Tensor) lying on the rotation axis or one of the following strings: 'bounds', 'individual'. Vector tensors are rotated about the origin if pivot is not given as a Tensor.

Returns

Rotated vector as Tensor

def rotation_angles(rot: phiml.math._tensors.Tensor)

Compute the scalar x in 2D or the Euler angles in 3D from a given rotation matrix. This function returns one valid solution but often, there are multiple solutions.

Args

rot: Rotation matrix as created by phi.math.rotation_matrix(). Must have exactly one channel and one dual dimension with equally-ordered elements.

Returns

Scalar x in 2D, Euler angles

def rotation_matrix(x: float | phiml.math._tensors.Tensor | None, matrix_dim=(vectorᶜ=None)) ‑> phiml.math._tensors.Tensor | None

Create a 2D or 3D rotation matrix from the corresponding angle(s).

Args

x:
2D: scalar angle
3D: Either vector pointing along the rotation axis with rotation angle as length or Euler angles.
Euler angles need to be laid out along a angle channel dimension with dimension names listing the spatial dimensions.
E.g. a 90° rotation about the z-axis is represented by vec('angles', x=0, y=0, z=PI/2).
If a rotation matrix is passed for angle, it is returned without modification.
matrix_dim: Matrix dimension for 2D rotations. In 3D, the channel dimension of angle is used.

Returns

Matrix containing matrix_dim in primal and dual form as well as all non-channel dimensions of x.

def rotation_matrix_from_axis_and_angle(axis: phiml.math._tensors.Tensor, angle: phiml.math._tensors.Tensor | float, vec_dim='vector', is_axis_normalized=False, epsilon=1e-05) ‑> phiml.math._tensors.Tensor

Computes a rotation matrix that rotates by angle around axis.

Args

axis: 3D vector. Tensor with channel dim called 'vector'.
angle: Rotation angle.
is_axis_normalized: Whether axis has length 1.
epsilon: Minimum axis length. For shorter axes, the unit matrix is returned.

Returns

Rotation matrix as Tensor with 'vector' dim and its dual counterpart.

def rotation_matrix_from_directions(source_dir: phiml.math._tensors.Tensor, target_dir: phiml.math._tensors.Tensor, vec_dim: str = 'vector', epsilon=None) ‑> phiml.math._tensors.Tensor

Computes a rotation matrix A, such that target_dir = A @ source_dir

Args

source_dir: Two or three-dimensional vector. Tensor with channel dim called 'vector'.
target_dir: Two or three-dimensional vector. Tensor with channel dim called 'vector'.

Returns

Rotation matrix as Tensor with 'vector' dim and its dual counterpart.

def sample_function(f: Callable, elements: phi.geom._geom.Geometry, at: str, extrapolation: phiml.math.extrapolation.Extrapolation) ‑> phiml.math._tensors.Tensor

Calls f, passing either the elements directly or the relevant sample points as a Tensor, depending on the signature of f.

Args

f: Function taking a Geometry or location Tensor´ and returning aTensor`. A Geometry will be passed if the first argument of f is called geometry or geo or ends with _geo.
elements: Geometry on which to sample f.
at: Set of sample points, see Geometry.sets.
extrapolation: Determines which boundary points are relevant.

Returns

Sampled values as Tensor.

def sample_sdf(geometry: phi.geom._geom.Geometry, bounds: phi.geom._box.BaseBox | phi.geom._grid.UniformGrid = None, resolution: phiml.math._shape.Shape = (), approximate_outside=False, rebuild: str | None = None, valid_dist=None, rel_margin=0.1, abs_margin=0.0, cache_surface=False, **resolution_: int) ‑> phi.geom._sdf_grid.SDFGrid

Build a grid of signed distance values for a given Geometry object.

Args

geometry: Geometry to capture.
bounds: Grid limits in world space.
resolution: Grid resolution.
**resolution_: Grid resolution as kwargs, e.g. x=64, y=32.
approximate_outside: Whether queries outside the SDF grid should return approximate values. This requires additional computations.
rebuild: If 'from-surface', SDF values are calculated from a narrow strip above the enclosed surface. This is more accurate but requires additional steps. If None (default), SDF values are queried from geometry. 'auto' rebuilds when geometry querying is expected to be in accurate.

Returns

SDF grid as Geometry.

def scale(obj: ~GeometricType, scale: phiml.math._tensors.Tensor | float, pivot: phiml.math._tensors.Tensor = None, dim='vector') ‑> ~GeometricType

Scale a Geometry or vector Tensor about a pivot point.

Args

obj: Geometry to scale.
scale: Scaling factor.
pivot: Point that stays fixed under the scaling operation. Defaults to the bounding box center.

Returns

Rotated Geometry

def stack(values: Sequence[~PhiTreeNodeType] | Dict[str, ~PhiTreeNodeType], dim: phiml.math._shape.Shape | str, expand_values=False, simplify=False, layout_non_matching=False, **kwargs) ‑> ~PhiTreeNodeType

Stacks values along the new dimension dim. All values must have the same spatial, instance and channel dimensions. If the dimension sizes vary, the resulting tensor will be non-uniform. Batch dimensions will be added as needed.

Stacking tensors is performed lazily, i.e. the memory is allocated only when needed. This makes repeated stacking and slicing along the same dimension very efficient, i.e. jit-compiled functions will not perform these operations.

Args

values: Collection of phiml.math.magic.Shapable, such as phiml.math.Tensor If a dict, keys must be of type str and are used as item names along dim.
dim: Shape with a least one dimension. None of these dimensions can be present with any of the values. If dim is a single-dimension shape, its size is determined from len(values) and can be left undefined (None). If dim is a multi-dimension shape, its volume must be equal to len(values).
expand_values: If True, will first add missing dimensions to all values, not just batch dimensions. This allows tensors with different dimensions to be stacked. The resulting tensor will have all dimensions that are present in values. If False, this may return a non-numeric object instead.
simplify: If True and all values are equal, returns one value without adding the dimension.
layout_non_matching: If non-matching values should be stacked using a Layout object, i.e. should be put into a named list instead.
**kwargs: Additional keyword arguments required by specific implementations. Adding spatial dimensions to fields requires the bounds: Box argument specifying the physical extent of the new dimensions. Adding batch dimensions must always work without keyword arguments.

Returns

Tensor containing values stacked along dim.

Examples

>>> stack({'x': 0, 'y': 1}, channel('vector'))
(x=0, y=1)

>>> stack([math.zeros(batch(b=2)), math.ones(batch(b=2))], channel(c='x,y'))
(x=0.000, y=1.000); (x=0.000, y=1.000) (bᵇ=2, cᶜ=x,y)

>>> stack([vec(x=1, y=0), vec(x=2, y=3.)], batch('b'))
(x=1.000, y=0.000); (x=2.000, y=3.000) (bᵇ=2, vectorᶜ=x,y)

def surface_mesh(geo: phi.geom._geom.Geometry, rel_dx: float = None, abs_dx: float = None, method='auto') ‑> phi.geom._mesh.Mesh

Create a surface Mesh from a Geometry.

Args

geo: Geometry to convert. Must implement approximate_signed_distance.
rel_dx: Relative mesh resolution as fraction of bounding box size.
abs_dx: Absolute mesh resolution. If both rel_dx and abs_dx are provided, the lower value is used.
method: 'auto' to select based on the type of geo. 'lewiner' or 'lorensen' for marching cubes.

Returns

Mesh if there is any geometry

def union(*geometries, dim=(unionⁱ=None)) ‑> phi.geom._geom.Geometry

Union of the given geometries. A point lies inside the union if it lies within at least one of the geometries.

Args

*geometries: arbitrary geometries with same spatial dims. Arbitrary batch dims are allowed.
dim: Union dimension. This must be an instance dimension.

Returns

union Geometry

Classes

class BaseBox

Abstract base type for box-like geometries.

Expand source code

class BaseBox(Geometry):  # not a Subwoofer
    """
    Abstract base type for box-like geometries.
    """

    def __ne__(self, other):
        return not self == other

    @property
    def shape(self):
        raise NotImplementedError()

    @property
    def center(self) -> Tensor:
        raise NotImplementedError()

    @property
    def size(self) -> Tensor:
        raise NotImplementedError(self)

    @property
    def half_size(self) -> Tensor:
        raise NotImplementedError(self)

    @property
    def lower(self) -> Tensor:
        raise NotImplementedError(self)

    @property
    def upper(self) -> Tensor:
        raise NotImplementedError(self)

    @property
    def rotation_matrix(self) -> Optional[Tensor]:
        raise NotImplementedError(self)

    @property
    def is_size_variable(self):
        raise NotImplementedError(self)

    @property
    def volume(self) -> Tensor:
        return math.prod(self.size, 'vector')

    def bounding_radius(self):
        return math.vec_length(self.half_size)

    def global_to_local(self, global_position: Tensor, scale=True, origin='lower') -> Tensor:
        """
        Transform world-space coordinates into box-space coordinates.

        Args:
            global_position: World-space coordinates.
            scale: Whether to re-scale the output so that [0, 1] or [-1, 1] represent the box for `origin='lower'` or `origin='center'`, respectively.
            origin: 'lower' or 'center'

        Returns:
            Box-space coordinate `Tensor`
        """
        assert origin in ['lower', 'center', 'upper']
        origin_loc = getattr(self, origin)
        pos = global_position if math.always_close(origin_loc, 0) else global_position - origin_loc
        pos = rotate(pos, self.rotation_matrix, invert=True)
        if scale:
            pos /= (self.half_size if origin == 'center' else self.size)
        return pos

    def local_to_global(self, local_position, scale=True, origin='lower'):
        assert origin in ['lower', 'center', 'upper']
        origin_loc = getattr(self, origin)
        pos = local_position * (self.half_size if origin == 'center' else self.size) if scale else local_position
        return rotate(pos, self.rotation_matrix) + origin_loc

    def largest(self, dim: DimFilter) -> 'BaseBox':
        dim = self.shape.without('vector').only(dim)
        if not dim:
            return self
        return Box(math.min(self.lower, dim), math.max(self.upper, dim))

    def smallest(self, dim: DimFilter) -> 'BaseBox':
        dim = self.shape.without('vector').only(dim)
        if not dim:
            return self
        return Box(math.max(self.lower, dim), math.min(self.upper, dim))

    def lies_inside(self, location: Tensor):
        assert self.rotation_matrix is None, f"Please override lies_inside() for rotated boxes"
        bool_inside = (location >= self.lower) & (location <= self.upper)
        bool_inside = math.all(bool_inside, 'vector')
        bool_inside = math.any(bool_inside, self.shape.instance - instance(location))  # union for instance dimensions
        return bool_inside

    def approximate_signed_distance(self, location: Union[Tensor, tuple]):
        """
        Computes the signed L-infinity norm (manhattan distance) from the location to the nearest side of the box.
        For an outside location `l` with the closest surface point `s`, the distance is `max(abs(l - s))`.
        For inside locations it is `-max(abs(l - s))`.

        Args:
          location: float tensor of shape (batch_size, ..., rank)

        Returns:
          float tensor of shape (*location.shape[:-1], 1).

        """
        location = self.global_to_local(location, scale=False, origin='center')
        distance = math.abs(location) - self.half_size
        distance = math.max(distance, 'vector')
        distance = math.min(distance, self.shape.instance)  # union for instance dimensions
        return distance

    def push(self, positions: Tensor, outward: bool = True, shift_amount: float = 0) -> Tensor:
        loc_to_center = self.global_to_local(positions, scale=False, origin='center')
        sgn_dist_from_surface = math.abs(loc_to_center) - self.half_size
        rotation_matrix = self.rotation_matrix
        if outward:
            # --- get negative distances (particles are inside) towards the nearest boundary and add shift_amount ---
            distances_of_interest = (sgn_dist_from_surface == math.max(sgn_dist_from_surface, 'vector')) & (sgn_dist_from_surface < 0)
            shift = distances_of_interest * (sgn_dist_from_surface - shift_amount)
            # ToDo reduce instance dim
        else:  # inward
            shift = (sgn_dist_from_surface + shift_amount) * (sgn_dist_from_surface > 0)  # get positive distances (particles are outside) and add shift_amount
            if instance(self):
                shift, loc_to_center, rotation_matrix = math.at_min((shift, loc_to_center, rotation_matrix), key=math.vec_length(shift), dim=instance)
            shift = math.where(abs(shift) > abs(loc_to_center), abs(loc_to_center), shift)  # ensure inward shift ends at center
        shift = rotate(shift, rotation_matrix)
        return positions + math.where(loc_to_center < 0, 1, -1) * shift

    def approximate_closest_surface(self, location: Tensor) -> Tuple[Tensor, Tensor, Tensor, Tensor, Tensor]:
        loc_to_center = self.global_to_local(location, scale=False, origin='center')
        sgn_surf_delta = math.abs(loc_to_center) - self.half_size
        if instance(self):
            raise NotImplementedError
            self_center, self_radius, sgn_dist, center_delta, center_dist = math.at_min((self.center, self.radius, sgn_dist, center_delta, center_dist), key=abs(sgn_dist), dim=instance)
        # is_inside = math.all(sgn_surf_delta < 0, 'vector')
        # abs_surf_delta = abs(sgn_surf_delta)
        max_sgn_dist = math.max(sgn_surf_delta, 'vector')
        normal_axis = max_sgn_dist == sgn_surf_delta  # ToDo only one if inside
        normal = math.vec_normalize(normal_axis * math.sign(loc_to_center))
        if self.rotation_matrix is not None:
            normal = rotate(normal, self.rotation_matrix)
        surf_to_center = math.where(normal_axis, math.sign(loc_to_center) * self.half_size, loc_to_center)
        closest_to_center = math.clip(surf_to_center, -self.half_size, self.half_size)
        surface_pos = self.local_to_global(closest_to_center, scale=False, origin='center')
        delta = surface_pos - location
        face_index = expand(0, non_channel(location))
        offset = normal.vector @ surface_pos.vector
        sgn_surf_dist = math.vec_length(delta) * math.sign(max_sgn_dist)
        return sgn_surf_dist, delta, normal, offset, face_index

    def project(self, *dimensions: str):
        """ Project this box into a lower-dimensional space. """
        warnings.warn("Box.project(dims) is deprecated. Use Box.vector[dims] instead", DeprecationWarning, stacklevel=2)
        return self.vector[dimensions]

    def sample_uniform(self, *shape: Shape) -> Tensor:
        uniform = math.random_uniform(self.shape.non_singleton.without('vector'), *shape, self.shape['vector'])
        return self.lower + uniform * self.size

    def sample_uniform_surface(self, *shape: Shape) -> Tensor:
        assert not instance(self), "sample_uniform_surface not yet supported for unions of boxes"
        samples = math.random_uniform(self.shape.non_singleton.non_instance, *shape, low=self.lower, high=self.upper)
        which = math.random_uniform(samples.shape.without('vector'))
        lo_or_up = math.where(which > .5, self.upper, self.lower)
        which = which * 2 % 1
        # --- which axis ---
        areas = self.face_areas
        total_area = math.sum(areas)
        frac_area = math.sum(areas / total_area, '~side')
        cum_area = math.cumulative_sum(frac_area, '~vector')
        axis = math.min(math.where(which <= cum_area, math.range(self.shape['vector'].as_dual()), self.spatial_rank), '~vector')
        axis_one_hot = math.scatter(math.zeros(samples.shape, dtype=bool), expand(axis, channel(index='vector')), True, treat_as_batch=samples.shape.without('vector'))
        math.assert_close(1, math.sum(axis_one_hot, 'vector'))
        samples = math.where(axis_one_hot, lo_or_up, samples)
        return samples

    def corner_representation(self) -> 'Box':
        assert self.rotation_matrix is None, f"corner_representation does not support rotations"
        return Box(self.lower, self.upper)

    box = corner_representation

    def center_representation(self, size_variable=True) -> 'Cuboid':
        return Cuboid(self.center, self.half_size, size_variable=size_variable)

    def contains(self, other: 'BaseBox'):
        """ Tests if the other box lies fully inside this box. """
        return np.all(other.lower >= self.lower) and np.all(other.upper <= self.upper)

    def scaled(self, factor: Union[float, Tensor]) -> 'BaseBox':
        return Cuboid(self.center, self.half_size * factor, size_variable=True)

    @property
    def boundary_elements(self) -> Dict[Any, Dict[str, slice]]:
        return {}

    @property
    def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
        return {}

    @property
    def faces(self) -> 'Geometry':
        return Cuboid(self.face_centers, self._half_size, self._rotation_matrix, size_variable=False)

    @property
    def face_centers(self) -> Tensor:
        return self.center + self.face_normals * self.half_size

    @property
    def face_normals(self) -> Tensor:
        unit_vectors = math.to_float(math.range(self.shape['vector']) == math.range(dual(**self.shape['vector'].untyped_dict)))
        vectors = rotate(unit_vectors, self.rotation_matrix)
        return vectors * math.vec(dual('side'), lower=-1, upper=1)

    @property
    def face_areas(self) -> Tensor:
        others_mask = math.range(self.shape['vector']) != math.range(dual(**self.shape['vector'].untyped_dict))
        result = math.exp(math.sum(math.log(self.size) * others_mask, 'vector'))
        return expand(result, dual(side='lower,upper'))  # ~vector

    @property
    def face_shape(self) -> Shape:
        return self.shape.without('vector') & dual(side='lower,upper') & dual(**self.shape['vector'].untyped_dict)

    @property
    def corners(self):
        to_face = self.face_normals[{'~side': 'upper'}] * math.rename_dims(self.half_size, 'vector', dual)
        lower_upper = math.meshgrid(math.dual, **{dim: [-1, 1] for dim in self.vector.item_names}, stack_dim=dual('vector'))  # (x=2, y=2, ... vector=x,y,...)
        to_corner = math.sum(lower_upper * to_face, '~vector')
        return self.center + to_corner

Ancestors

phi.geom._geom.Geometry

Subclasses

phi.geom._box.Box
phi.geom._box.Cuboid
phi.geom._grid.UniformGrid

Instance variables

prop boundary_elements : Dict[Any, Dict[str, slice]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[Any, Dict[str, slice]]:
    return {}

prop boundary_faces : Dict[Any, Dict[str, slice]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
    return {}

prop center : phiml.math._tensors.Tensor

Center location in single channel dimension.

Expand source code

@property
def center(self) -> Tensor:
    raise NotImplementedError()

prop corners

Returns

Corner locations as phiml.math.Tensor. Corners belonging to one object or cell are listed along dual dimensions. If the object has no corners, a size-0 tensor with the correct vector and instance dims is returned.

Expand source code

@property
def corners(self):
    to_face = self.face_normals[{'~side': 'upper'}] * math.rename_dims(self.half_size, 'vector', dual)
    lower_upper = math.meshgrid(math.dual, **{dim: [-1, 1] for dim in self.vector.item_names}, stack_dim=dual('vector'))  # (x=2, y=2, ... vector=x,y,...)
    to_corner = math.sum(lower_upper * to_face, '~vector')
    return self.center + to_corner

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    others_mask = math.range(self.shape['vector']) != math.range(dual(**self.shape['vector'].untyped_dict))
    result = math.exp(math.sum(math.log(self.size) * others_mask, 'vector'))
    return expand(result, dual(side='lower,upper'))  # ~vector

prop face_centers : phiml.math._tensors.Tensor

Center of face connecting a pair of cells. Shape (elements, ~, vector). Here, ~ represents arbitrary internal dual dimensions, such as ~staggered_direction or ~elements. Returns 0-vectors for unconnected cells.

Expand source code

@property
def face_centers(self) -> Tensor:
    return self.center + self.face_normals * self.half_size

prop face_normals : phiml.math._tensors.Tensor

Normal vectors of cell faces, including boundary faces. Shape (elements, ~, vector). For meshes, The vectors point out of the primal cells and into the dual cells.

Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.

Expand source code

@property
def face_normals(self) -> Tensor:
    unit_vectors = math.to_float(math.range(self.shape['vector']) == math.range(dual(**self.shape['vector'].untyped_dict)))
    vectors = rotate(unit_vectors, self.rotation_matrix)
    return vectors * math.vec(dual('side'), lower=-1, upper=1)

prop face_shape : phiml.math._shape.Shape

Returns

Full Shape to identify each face of this Geometry, including instance/spatial dimensions for the elements and dual dimensions listing the faces per element. If this Geometry has no faces, returns an empty Shape.

Expand source code

@property
def face_shape(self) -> Shape:
    return self.shape.without('vector') & dual(side='lower,upper') & dual(**self.shape['vector'].untyped_dict)

prop faces : Geometry

Expand source code

@property
def faces(self) -> 'Geometry':
    return Cuboid(self.face_centers, self._half_size, self._rotation_matrix, size_variable=False)

prop half_size : phiml.math._tensors.Tensor

Expand source code

@property
def half_size(self) -> Tensor:
    raise NotImplementedError(self)

prop is_size_variable

Expand source code

@property
def is_size_variable(self):
    raise NotImplementedError(self)

prop lower : phiml.math._tensors.Tensor

Expand source code

@property
def lower(self) -> Tensor:
    raise NotImplementedError(self)

prop rotation_matrix : phiml.math._tensors.Tensor | None

Expand source code

@property
def rotation_matrix(self) -> Optional[Tensor]:
    raise NotImplementedError(self)

prop shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self):
    raise NotImplementedError()

prop size : phiml.math._tensors.Tensor

Expand source code

@property
def size(self) -> Tensor:
    raise NotImplementedError(self)

prop upper : phiml.math._tensors.Tensor

Expand source code

@property
def upper(self) -> Tensor:
    raise NotImplementedError(self)

prop volume : phiml.math._tensors.Tensor

phi.math.Tensor representing the volume of each element. The result retains batch, spatial and instance dimensions.

Expand source code

@property
def volume(self) -> Tensor:
    return math.prod(self.size, 'vector')

Methods

def approximate_closest_surface(self, location: phiml.math._tensors.Tensor) ‑> Tuple[phiml.math._tensors.Tensor, phiml.math._tensors.Tensor, phiml.math._tensors.Tensor, phiml.math._tensors.Tensor, phiml.math._tensors.Tensor]

Find the closest surface face of this geometry given a point that can be outside or inside the geometry.

Args

location: Tensor with a single channel dimension called vector. Can have arbitrary other dimensions.

Returns

signed_distance: Scalar signed distance from location to the closest point on the surface. Positive values indicate the point lies outside the geometry, negative values indicate the point lies inside the geometry.
delta: Vector-valued distance vector from location to the closest point on the surface.
normal: Closest surface normal vector.
offset: Min distance of a surface-tangential plane from 0 as a scalar.
face_index: (Optional) An index vector pointing at the closest face.

def approximate_signed_distance(self, location: phiml.math._tensors.Tensor | tuple)

Computes the signed L-infinity norm (manhattan distance) from the location to the nearest side of the box. For an outside location l with the closest surface point s, the distance is max(abs(l - s)). For inside locations it is -max(abs(l - s)).

Args

location: float tensor of shape (batch_size, …, rank)

Returns

float tensor of shape (*location.shape[:-1], 1).

def bounding_radius(self)

Returns the radius of a Sphere object that fully encloses this geometry. The sphere is centered at the center of this geometry.

If this geometry consists of multiple parts listed along instance/spatial dims, these dims are retained, giving the bounds of each part. If these dims are not present on the result, all parts are assumed to have the same bounds.

def box(self) ‑> phi.geom._box.Box

def center_representation(self, size_variable=True) ‑> phi.geom._box.Cuboid

def contains(self, other: BaseBox)

Tests if the other box lies fully inside this box.

def corner_representation(self) ‑> phi.geom._box.Box

def global_to_local(self, global_position: phiml.math._tensors.Tensor, scale=True, origin='lower') ‑> phiml.math._tensors.Tensor

Transform world-space coordinates into box-space coordinates.

Args

global_position: World-space coordinates.
scale: Whether to re-scale the output so that [0, 1] or [-1, 1] represent the box for origin='lower' or origin='center', respectively.
origin: 'lower' or 'center'

Returns

Box-space coordinate Tensor

def lies_inside(self, location: phiml.math._tensors.Tensor)

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def local_to_global(self, local_position, scale=True, origin='lower')

def project(self, *dimensions: str)

Project this box into a lower-dimensional space.

def push(self, positions: phiml.math._tensors.Tensor, outward: bool = True, shift_amount: float = 0) ‑> phiml.math._tensors.Tensor

Shifts positions either into or out of geometry.

Args

positions: Tensor holding positions to shift
outward: Flag for indicating inward (False) or outward (True) shift
shift_amount: Minimum distance between positions and surface after shifting.

Returns

Tensor holding shifted positions.

def sample_uniform(self, *shape: phiml.math._shape.Shape) ‑> phiml.math._tensors.Tensor

Samples uniformly distributed random points inside this volume.

Args

*shape: How many points to sample per individual geometry.

Returns

Tensor containing all dimensions from Geometry.shape, shape as well as a channel dimension vector matching the dimensionality of this Geometry.

def sample_uniform_surface(self, *shape: phiml.math._shape.Shape) ‑> phiml.math._tensors.Tensor

def scaled(self, factor: phiml.math._tensors.Tensor | float) ‑> phi.geom._box.BaseBox

Scales each individual geometry by factor. The individual center points act as pivots for the operation.

Args

factor: Returns:

class Box (lower: phiml.math._tensors.Tensor = None, upper: phiml.math._tensors.Tensor = None, **size: float | phiml.math._tensors.Tensor | tuple | list | None)

Simple cuboid defined by location of lower and upper corner in physical space.

Boxes can be constructed either from two positional vector arguments (lower, upper) or by specifying the limits by dimension name as kwargs.

Examples

>>> Box(x=1, y=1)  # creates a two-dimensional unit box with `lower=(0, 0)` and `upper=(1, 1)`.
>>> Box(x=(None, 1), y=(0, None)  # creates a Box with `lower=(-inf, 0)` and `upper=(1, inf)`.

The slicing constructor was updated in version 2.2 and now requires the dimension order as the first argument.

>>> Box['x,y', 0:1, 0:1]  # creates a two-dimensional unit box with `lower=(0, 0)` and `upper=(1, 1)`.
>>> Box['x,y', :1, 0:]  # creates a Box with `lower=(-inf, 0)` and `upper=(1, inf)`.

Args

lower: physical location of lower corner
upper: physical location of upper corner
**size: Specify size by dimension, either as int or tuple containing (lower, upper).

Expand source code

class Box(BaseBox, metaclass=BoxType):
    """
    Simple cuboid defined by location of lower and upper corner in physical space.

    Boxes can be constructed either from two positional vector arguments `(lower, upper)` or by specifying the limits by dimension name as `kwargs`.

    Examples:
        >>> Box(x=1, y=1)  # creates a two-dimensional unit box with `lower=(0, 0)` and `upper=(1, 1)`.
        >>> Box(x=(None, 1), y=(0, None)  # creates a Box with `lower=(-inf, 0)` and `upper=(1, inf)`.

        The slicing constructor was updated in version 2.2 and now requires the dimension order as the first argument.

        >>> Box['x,y', 0:1, 0:1]  # creates a two-dimensional unit box with `lower=(0, 0)` and `upper=(1, 1)`.
        >>> Box['x,y', :1, 0:]  # creates a Box with `lower=(-inf, 0)` and `upper=(1, inf)`.
    """

    def __init__(self, lower: Tensor = None, upper: Tensor = None, **size: Optional[Union[float, Tensor, tuple, list]]):
        """
        Args:
          lower: physical location of lower corner
          upper: physical location of upper corner
          **size: Specify size by dimension, either as `int` or `tuple` containing (lower, upper).
        """
        if lower is not None:
            assert isinstance(lower, Tensor), f"lower must be a Tensor but got {type(lower)}"
            assert 'vector' in lower.shape, "lower must have a vector dimension"
            assert lower.vector.item_names is not None, "vector dimension of lower must list spatial dimension order"
            self._lower = lower
        if upper is not None:
            assert isinstance(upper, Tensor), f"upper must be a Tensor but got {type(upper)}"
            assert 'vector' in upper.shape, "lower must have a vector dimension"
            assert upper.vector.item_names is not None, "vector dimension of lower must list spatial dimension order"
            self._upper = upper
        else:
            lower = []
            upper = []
            for item in size.values():
                if isinstance(item, (tuple, list)):
                    assert len(item) == 2, f"Box kwargs must be either dim=upper or dim=(lower,upper) but got {item}"
                    lo, up = item
                    lower.append(lo)
                    upper.append(up)
                elif item is None:
                    lower.append(-INF)
                    upper.append(INF)
                else:
                    lower.append(0)
                    upper.append(item)
            lower = [-INF if l is None else l for l in lower]
            upper = [INF if u is None else u for u in upper]
            self._upper = math.wrap(upper, math.channel(vector=tuple(size.keys())))
            self._lower = math.wrap(lower, math.channel(vector=tuple(size.keys())))
        vector_shape = self._lower.shape & self._upper.shape
        self._lower = math.expand(self._lower, vector_shape)
        self._upper = math.expand(self._upper, vector_shape)
        if self.size.vector.item_names is None:
            warnings.warn("Creating a Box without item names prevents certain operations like project()", DeprecationWarning, stacklevel=2)
        self._shape = self._lower.shape & self._upper.shape

    def __getitem__(self, item):
        item = _keep_vector(slicing_dict(self, item))
        return Box(self._lower[item], self._upper[item])

    @staticmethod
    def __stack__(values: tuple, dim: Shape, **kwargs) -> 'Geometry':
        if all(isinstance(v, Box) for v in values):
            return NotImplemented  # stack attributes
        else:
            return Geometry.__stack__(values, dim, **kwargs)

    def without(self, dims: Tuple[str, ...]):
        remaining = list(self.shape.get_item_names('vector'))
        for dim in dims:
            if dim in remaining:
                remaining.remove(dim)
        return self.vector[remaining]

    def __variable_attrs__(self):
        return '_lower', '_upper'

    def __value_attrs__(self):
        return '_lower', '_upper'

    @property
    def shape(self):
        if self._lower is None or self._upper is None:
            return self._shape
        return self._lower.shape & self._upper.shape

    @property
    def lower(self):
        return self._lower

    @property
    def upper(self):
        return self._upper

    @property
    def size(self):
        return self.upper - self.lower

    @property
    def center(self):
        return 0.5 * (self.lower + self.upper)

    @property
    def half_size(self):
        return self.size * 0.5

    @property
    def rotation_matrix(self) -> Optional[Tensor]:
        return None

    @property
    def is_size_variable(self):
        raise False

    def at(self, center: Tensor) -> 'BaseBox':
        return Cuboid(center, self.half_size, self.rotation_matrix)

    def shifted(self, delta, **delta_by_dim) -> 'Box':
        return Box(self.lower + delta, self.upper + delta)

    def rotated(self, angle) -> 'Cuboid':
        return self.center_representation().rotated(angle)

    def __mul__(self, other):
        if not isinstance(other, Box):
            return NotImplemented
        lower = self._lower.vector.unstack(self.spatial_rank) + other._lower.vector.unstack(other.spatial_rank)
        upper = self._upper.vector.unstack(self.spatial_rank) + other._upper.vector.unstack(other.spatial_rank)
        names = self._upper.vector.item_names + other._upper.vector.item_names
        lower = math.stack(lower, math.channel(vector=names))
        upper = math.stack(upper, math.channel(vector=names))
        return Box(lower, upper)

    def bounding_half_extent(self):
        return self.half_size

    def __repr__(self):
        if self._lower is None or self._upper is None:  # traced
            return f"Box[traced, shape={self._shape}]"
        if self.shape.non_channel.volume == 1:
            item_names = self.size.vector.item_names
            if item_names:
                return f"Box({', '.join([f'{dim}=({lo}, {up})' for dim, lo, up in zip(item_names, self._lower, self._upper)])})"
            else:  # deprecated
                return 'Box[%s at %s]' % ('x'.join([str(x) for x in self.size.numpy().flatten()]), ','.join([str(x) for x in self.lower.numpy().flatten()]))
        else:
            return f'Box[shape={self.shape}]'

Ancestors

phi.geom._box.BaseBox
phi.geom._geom.Geometry

Instance variables

prop center

Center location in single channel dimension.

Expand source code

@property
def center(self):
    return 0.5 * (self.lower + self.upper)

prop half_size

Expand source code

@property
def half_size(self):
    return self.size * 0.5

prop is_size_variable

Expand source code

@property
def is_size_variable(self):
    raise False

prop lower

Expand source code

@property
def lower(self):
    return self._lower

prop rotation_matrix : phiml.math._tensors.Tensor | None

Expand source code

@property
def rotation_matrix(self) -> Optional[Tensor]:
    return None

prop shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self):
    if self._lower is None or self._upper is None:
        return self._shape
    return self._lower.shape & self._upper.shape

prop size

Expand source code

@property
def size(self):
    return self.upper - self.lower

prop upper

Expand source code

@property
def upper(self):
    return self._upper

Methods

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._box.BaseBox

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_half_extent(self)

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

Let the bounding half-extent have value e in dimension d (extent[...,d] = e). Then, no point of the geometry lies further away from its center point than e along d (in both axis directions).

def rotated(self, angle) ‑> phi.geom._box.Cuboid

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def shifted(self, delta, **delta_by_dim) ‑> phi.geom._box.Box

Returns a translated version of this geometry.

See Also: Geometry.at().

Args

delta: direction vector
delta: Tensor:

Returns

Geometry: shifted geometry

def without(self, dims: Tuple[str, ...])

class Cuboid (center: phiml.math._tensors.Tensor = 0, half_size: phiml.math._tensors.Tensor | float = None, rotation: phiml.math._tensors.Tensor | None = None, size_variable=True, **size: phiml.math._tensors.Tensor | float)

Box specified by center position and half size.

Args

center: Center position
half_size: Half-size of the cuboid as vector or scalar
rotation: Rotation angle(s) or rotation matrix.
**size: Alternative way of specifying the size. If used, half_size must not be specified.

Expand source code

class Cuboid(BaseBox):
    """Box specified by center position and half size."""

    def __init__(self,
                 center: Tensor = 0,
                 half_size: Union[float, Tensor] = None,
                 rotation: Optional[Tensor] = None,
                 size_variable=True,
                 **size: Union[float, Tensor]):
        """
        Args:
            center: Center position
            half_size: Half-size of the cuboid as vector or scalar
            rotation: Rotation angle(s) or rotation matrix.
            **size: Alternative way of specifying the size. If used, `half_size` must not be specified.
        """
        if half_size is not None:
            assert isinstance(half_size, Tensor), "half_size must be a Tensor"
            assert 'vector' in half_size.shape, f"Cuboid size must have a 'vector' dimension."
            assert half_size.shape.get_item_names('vector') is not None, f"Vector dimension must list spatial dimensions as item names. Use the syntax Cuboid(x=x, y=y) to assign names."
            self._half_size = half_size
        else:
            self._half_size = math.wrap(tuple(size.values()), math.channel(vector=tuple(size.keys()))) * 0.5
        center = wrap(center)
        if 'vector' not in center.shape or center.shape.get_item_names('vector') is None:
            center = math.expand(center, channel(self._half_size))
        self._center = center
        self._rotation_matrix = None if rotation is None else rotation_matrix(rotation, matrix_dim=center.shape['vector'])
        self._size_variable = size_variable

    def __repr__(self):
        return f"Cuboid(center={self._center}, half_size={self._half_size})"

    def __getitem__(self, item) -> 'Cuboid':
        item = _keep_vector(slicing_dict(self, item))
        return Cuboid(self._center[item], self._half_size[item], math.slice(self._rotation_matrix, item), size_variable=self._size_variable)

    @staticmethod
    def __stack__(values: tuple, dim: Shape, **kwargs) -> 'Geometry':
        if all(isinstance(v, Cuboid) for v in values):
            size_variable = any([c._size_variable for c in values])
            if any(v._rotation_matrix is not None for v in values):
                matrices = [v._rotation_matrix for v in values]
                if any(m is None for m in matrices):
                    any_angle = math.rotation_angles([m for m in matrices if m is not None][0])
                    unit_matrix = math.rotation_matrix(any_angle * 0)
                    matrices = [unit_matrix if m is None else m for m in matrices]
                rotation = stack(matrices, dim, **kwargs)
            else:
                rotation = None
            return Cuboid(stack([v.center for v in values], dim, **kwargs), stack([v.half_size for v in values], dim, **kwargs), rotation, size_variable=size_variable)
        else:
            return Geometry.__stack__(values, dim, **kwargs)

    def __variable_attrs__(self):
        return ('_center', '_half_size', '_rotation_matrix') if self._size_variable else ('_center', '_rotation_matrix')

    def __value_attrs__(self):
        return '_center',

    @property
    def center(self):
        return self._center

    @property
    def half_size(self):
        return self._half_size

    @property
    def shape(self):
        if self._center is None or self._half_size is None:
            return None
        return self._center.shape & self._half_size.shape

    @property
    def size(self):
        return 2 * self._half_size

    @property
    def lower(self):
        return self._center - self._half_size

    @property
    def upper(self):
        return self._center + self._half_size

    @property
    def rotation_matrix(self) -> Optional[Tensor]:
        return self._rotation_matrix

    @property
    def is_size_variable(self):
        return self._size_variable

    def at(self, center: Tensor) -> 'Cuboid':
        return Cuboid(center, self.half_size, self.rotation_matrix, size_variable=self._size_variable)

    def rotated(self, angle) -> 'Cuboid':
        if self._rotation_matrix is None:
            return Cuboid(self._center, self._half_size, angle, size_variable=self._size_variable)
        else:
            matrix = self._rotation_matrix @ (angle if dual(angle) else math.rotation_matrix(angle))
            return Cuboid(self._center, self._half_size, matrix, size_variable=self._size_variable)

    def bounding_half_extent(self) -> Tensor:
        if self._rotation_matrix is not None:
            to_face = self.face_normals[{'~side': 0}] * math.rename_dims(self._half_size, 'vector', dual)
            return math.sum(abs(to_face), '~vector')
        return self.half_size

    def lies_inside(self, location: Tensor) -> Tensor:
        location = self.global_to_local(location, scale=False, origin='center')  # scale can only be performed for finite sizes
        bool_inside = abs(location) <= self._half_size
        bool_inside = math.all(bool_inside, 'vector')
        bool_inside = math.any(bool_inside, self.shape.instance - instance(location))  # union for instance dimensions
        return bool_inside

Ancestors

phi.geom._box.BaseBox
phi.geom._geom.Geometry

Instance variables

prop center

Center location in single channel dimension.

Expand source code

@property
def center(self):
    return self._center

prop half_size

Expand source code

@property
def half_size(self):
    return self._half_size

prop is_size_variable

Expand source code

@property
def is_size_variable(self):
    return self._size_variable

prop lower

Expand source code

@property
def lower(self):
    return self._center - self._half_size

prop rotation_matrix : phiml.math._tensors.Tensor | None

Expand source code

@property
def rotation_matrix(self) -> Optional[Tensor]:
    return self._rotation_matrix

prop shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self):
    if self._center is None or self._half_size is None:
        return None
    return self._center.shape & self._half_size.shape

prop size

Expand source code

@property
def size(self):
    return 2 * self._half_size

prop upper

Expand source code

@property
def upper(self):
    return self._center + self._half_size

Methods

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._box.Cuboid

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_half_extent(self) ‑> phiml.math._tensors.Tensor

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def lies_inside(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def rotated(self, angle) ‑> phi.geom._box.Cuboid

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

class Cylinder (_center: phiml.math._tensors.Tensor, radius: phiml.math._tensors.Tensor, depth: phiml.math._tensors.Tensor, rotation: phiml.math._tensors.Tensor, axis: str, variable_attrs: Tuple[str, ...] = ('_center', 'radius', 'depth', 'rotation'), value_attrs: Tuple[str, ...] = ())

N-dimensional cylinder. Defined by center position, radius, depth, alignment axis, rotation.

For cylinders whose bottom and top lie outside the domain or are otherwise not needed, you may use infinite_cylinder() instead, which simplifies computations.

Expand source code

@dataclass(frozen=True)
class Cylinder(Geometry):
    """
    N-dimensional cylinder.
    Defined by center position, radius, depth, alignment axis, rotation.

    For cylinders whose bottom and top lie outside the domain or are otherwise not needed, you may use `infinite_cylinder` instead, which simplifies computations.
    """

    _center: Tensor
    radius: Tensor
    depth: Tensor
    rotation: Tensor  # rotation matrix
    axis: str

    variable_attrs: Tuple[str, ...] = ('_center', 'radius', 'depth', 'rotation')
    value_attrs: Tuple[str, ...] = ()

    @property
    def center(self) -> Tensor:
        return self._center

    @cached_property
    def shape(self) -> Shape:
        return self._center.shape & self.radius.shape & self.depth.shape & batch(self.rotation)

    @cached_property
    def radial_axes(self) -> Sequence[str]:
        return [d for d in self._center.vector.item_names if d != self.axis]

    @cached_property
    def volume(self) -> math.Tensor:
        return Sphere.volume_from_radius(self.radius, self.spatial_rank - 1) * self.depth

    @cached_property
    def up(self):
        return rotate(vec(**{d: 1 if d == self.axis else 0 for d in self._center.vector.item_names}), self.rotation)

    def with_radius(self, radius: Tensor) -> 'Cylinder':
        return Cylinder(self._center, wrap(radius), self.depth, self.rotation, self.axis, self.variable_attrs, self.value_attrs)

    def with_depth(self, depth: Tensor) -> 'Cylinder':
        return Cylinder(self._center, self.radius, wrap(depth), self.rotation, self.axis, self.variable_attrs, self.value_attrs)

    def lies_inside(self, location):
        pos = rotate(location - self._center, self.rotation, invert=True)
        r = pos.vector[self.radial_axes]
        h = pos.vector[self.axis]
        inside = (vec_squared(r) <= self.radius**2) & (h >= -.5*self.depth) & (h <= .5*self.depth)
        return math.any(inside, instance(self))  # union for instance dimensions

    def approximate_signed_distance(self, location: Union[Tensor, tuple]):
        location = rotate(location - self._center, self.rotation, invert=True)
        r = location.vector[self.radial_axes]
        h = location.vector[self.axis]
        top_h = .5*self.depth
        # bot_h = -.5*self.depth
        # --- Compute distances ---
        radial_outward = normalize(r, 'vector', epsilon=1e-5)
        surf_r = radial_outward * self.radius
        radial_dist2 = vec_squared(r)
        inside_cyl = radial_dist2 <= self.radius**2
        clamped_r = where(inside_cyl, r, surf_r)
        # --- Closest point on bottom / top ---
        sgn_dist_side = abs(h) - top_h
        # --- Closest point on cylinder ---
        sgn_dist_cyl = length(r, 'vector') - self.radius
        # inside (all <= 0) -> largest SDF, outside (any > 0) -> largest positive SDF
        sgn_dist = maximum(sgn_dist_cyl, sgn_dist_side)
        return math.min(sgn_dist, instance(self))

    def approximate_closest_surface(self, location: Tensor):
        location = rotate(location - self._center, self.rotation, invert=True)
        r = location.vector[self.radial_axes]
        h = location.vector[self.axis]
        top_h = .5*self.depth
        bot_h = -.5*self.depth
        # --- Compute distances ---
        radial_outward = normalize(r, 'vector', epsilon=1e-5)
        surf_r = radial_outward * self.radius
        radial_dist2 = vec_squared(r)
        inside_cyl = radial_dist2 <= self.radius**2
        clamped_r = where(inside_cyl, r, surf_r)
        # --- Closest point on bottom / top ---
        above = h >= 0
        flat_h = where(above, top_h, bot_h)
        on_flat = ncat([flat_h, clamped_r], self._center.shape['vector'])
        normal_flat = where(above, self.up, -self.up)
        # --- Closest point on cylinder ---
        clamped_h = clip(h, bot_h, top_h)
        on_cyl = ncat([surf_r, clamped_h], self._center.shape['vector'])
        normal_cyl = ncat([radial_outward, 0], self._center.shape['vector'], expand_values=True)
        # --- Choose closest ---
        d_flat = length(on_flat - location, 'vector')
        d_cyl = length(on_cyl - location, 'vector')
        flat_closer = d_flat <= d_cyl
        surf_point = where(flat_closer, on_flat, on_cyl)
        inside = inside_cyl & (h >= bot_h) & (h <= top_h)
        sgn_dist = minimum(d_flat, d_cyl) * where(inside, -1, 1)
        delta = surf_point - location
        normal = where(flat_closer, normal_flat, normal_cyl)
        delta = rotate(delta, self.rotation)
        normal = rotate(normal, self.rotation)
        idx = None
        if instance(self):
            sgn_dist, delta, normal, idx = math.min((sgn_dist, delta, normal, range), instance(self), key=sgn_dist)
        return sgn_dist, delta, normal, None, idx

    def sample_uniform(self, *shape: math.Shape):
        r = Sphere(self._center[self.radial_axes], self.radius).sample_uniform(*shape)
        h = math.random_uniform(*shape, -.5*self.depth, .5*self.depth)
        rh = ncat([r, h], self._center.shape['vector'])
        return rotate(rh, self.rotation)

    def bounding_radius(self):
        return length(vec(rad=self.radius, dep=.5*self.depth), 'vector')

    def bounding_half_extent(self, epsilon=1e-5):
        if self.rotation is not None:
            tip = abs(self.up) * .5 * self.depth
            return tip + self.radius * sqrt(maximum(epsilon, 1 - self.up**2))
        return ncat([.5*self.depth, expand(self.radius, channel(vector=self.radial_axes))], self._center.shape['vector'], expand_values=True)

    def at(self, center: Tensor) -> 'Geometry':
        return Cylinder(center, self.radius, self.depth, self.rotation, self.axis, self.variable_attrs, self.value_attrs)

    def rotated(self, angle):
        rot = self.rotation @ rotation_matrix(angle) if self.rotation is not None else rotation_matrix(angle)
        return Cylinder(self.center, self.radius, self.depth, rot, self.axis, self.variable_attrs, self.value_attrs)

    def scaled(self, factor: Union[float, Tensor]) -> 'Geometry':
        return Cylinder(self._center, self.radius * factor, self.depth * factor, self.rotation, self.axis, self.variable_attrs, self.value_attrs)

    def __getitem__(self, item):
        return getitem_dataclass(self, item, keepdims='vector')

    @staticmethod
    def __stack__(values: tuple, dim: Shape, **kwargs) -> 'Geometry':
        if all(isinstance(v, Cylinder) for v in values) and all(v.axis == values[0].axis for v in values):
            var_attrs = set()
            var_attrs.update(*[set(v.variable_attrs) for v in values])
            val_attrs = set()
            val_attrs.update(*[set(v.value_attrs) for v in values])
            if any(v.rotation is not None for v in values):
                matrices = [v.rotation for v in values]
                if any(m is None for m in matrices):
                    any_angle = rotation_angles([m for m in matrices if m is not None][0])
                    unit_matrix = rotation_matrix(any_angle * 0)
                    matrices = [unit_matrix if m is None else m for m in matrices]
                rotation = stack(matrices, dim, **kwargs)
            else:
                rotation = None
            center = stack([v.center for v in values], dim, simplify=True, **kwargs)
            radius = stack([v.radius for v in values], dim, simplify=True, **kwargs)
            depth = stack([v.depth for v in values], dim, simplify=True, **kwargs)
            return Cylinder(center, radius, depth, rotation, values[0].axis, tuple(var_attrs), tuple(val_attrs))
        else:
            return Geometry.__stack__(values, dim, **kwargs)

    @property
    def faces(self) -> 'Geometry':
        raise NotImplementedError(f"Cylinder.faces not implemented.")

    @property
    def face_centers(self) -> Tensor:
        raise NotImplementedError

    @property
    def face_areas(self) -> Tensor:
        flat = Sphere.volume_from_radius(self.radius, self.spatial_rank - 1)
        lateral = 2*PI*self.radius * self.depth
        return stack({'bottom': flat, 'top': flat, 'lateral': lateral}, dual('shell'), expand_values=True)

    @property
    def face_normals(self) -> Tensor:
        raise NotImplementedError

    @property
    def boundary_elements(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
        return {}

    @property
    def boundary_faces(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
        return {}

    @property
    def face_shape(self) -> Shape:
        return self.shape.without('vector') & dual(shell='bottom,top,lateral')

    @property
    def corners(self) -> Tensor:
        return math.zeros(self.shape & dual(corners=0))

    def __eq__(self, other):
        return Geometry.__eq__(self, other)

    def vertex_rings(self, count: Shape) -> Tensor:
        if self.spatial_rank == 3:
            angle = linspace(0, 2*PI, count)
            h = stack({'bot': -.5 * self.depth, 'top': .5 * self.depth}, '~face')
            s = sin(angle) * self.radius
            c = cos(angle) * self.radius
            r = stack([s, c], channel(vector=self.radial_axes))
            x = ncat([h, r], self._center.shape['vector'], expand_values=True)
            return rotate(x, self.rotation) + self._center
        raise NotImplementedError

Ancestors

phi.geom._geom.Geometry

Class variables

var axis : str
var depth : phiml.math._tensors.Tensor
var radius : phiml.math._tensors.Tensor
var rotation : phiml.math._tensors.Tensor
var value_attrs : Tuple[str, ...]
var variable_attrs : Tuple[str, ...]

Instance variables

prop boundary_elements : Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
    return {}

prop boundary_faces : Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
    return {}

prop center : phiml.math._tensors.Tensor

Center location in single channel dimension.

Expand source code

@property
def center(self) -> Tensor:
    return self._center

prop corners : phiml.math._tensors.Tensor

Returns

Expand source code

@property
def corners(self) -> Tensor:
    return math.zeros(self.shape & dual(corners=0))

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    flat = Sphere.volume_from_radius(self.radius, self.spatial_rank - 1)
    lateral = 2*PI*self.radius * self.depth
    return stack({'bottom': flat, 'top': flat, 'lateral': lateral}, dual('shell'), expand_values=True)

prop face_centers : phiml.math._tensors.Tensor

Expand source code

@property
def face_centers(self) -> Tensor:
    raise NotImplementedError

prop face_normals : phiml.math._tensors.Tensor

Normal vectors of cell faces, including boundary faces. Shape (elements, ~, vector). For meshes, The vectors point out of the primal cells and into the dual cells.

Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.

Expand source code

@property
def face_normals(self) -> Tensor:
    raise NotImplementedError

prop face_shape : phiml.math._shape.Shape

Returns

Expand source code

@property
def face_shape(self) -> Shape:
    return self.shape.without('vector') & dual(shell='bottom,top,lateral')

prop faces : Geometry

Expand source code

@property
def faces(self) -> 'Geometry':
    raise NotImplementedError(f"Cylinder.faces not implemented.")

var radial_axes

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

var shape

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

var up

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

var volume

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

Methods

def approximate_closest_surface(self, location: phiml.math._tensors.Tensor)

Find the closest surface face of this geometry given a point that can be outside or inside the geometry.

Args

location: Tensor with a single channel dimension called vector. Can have arbitrary other dimensions.

Returns

signed_distance: Scalar signed distance from location to the closest point on the surface. Positive values indicate the point lies outside the geometry, negative values indicate the point lies inside the geometry.
delta: Vector-valued distance vector from location to the closest point on the surface.
normal: Closest surface normal vector.
offset: Min distance of a surface-tangential plane from 0 as a scalar.
face_index: (Optional) An index vector pointing at the closest face.

def approximate_signed_distance(self, location: phiml.math._tensors.Tensor | tuple)

Computes the approximate distance from location to the surface of the geometry. Locations outside return positive values, inside negative values and zero exactly at the boundary.

The exact distance metric used depends on the geometry. The approximation holds close to the surface and the distance grows to infinity as the location is moved infinitely far from the geometry. The distance metric is differentiable and its gradients are bounded at every point in space.

When dealing with unions or collections of geometries (instance dimensions), the shortest distance to any instance is returned. This also holds for negative distances.

Args

location: Tensor with one channel dim vector matching the geometry's vector dim.

Returns

Float Tensor

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_half_extent(self, epsilon=1e-05)

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def bounding_radius(self)

Returns the radius of a Sphere object that fully encloses this geometry. The sphere is centered at the center of this geometry.

def lies_inside(self, location)

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def rotated(self, angle)

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def sample_uniform(self, *shape: phiml.math._shape.Shape)

Samples uniformly distributed random points inside this volume.

Args

*shape: How many points to sample per individual geometry.

Returns

Tensor containing all dimensions from Geometry.shape, shape as well as a channel dimension vector matching the dimensionality of this Geometry.

def scaled(self, factor: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Scales each individual geometry by factor. The individual center points act as pivots for the operation.

Args

factor: Returns:

def vertex_rings(self, count: phiml.math._shape.Shape) ‑> phiml.math._tensors.Tensor

def with_depth(self, depth: phiml.math._tensors.Tensor) ‑> phi.geom._cylinder.Cylinder

def with_radius(self, radius: phiml.math._tensors.Tensor) ‑> phi.geom._cylinder.Cylinder

class Geometry

Abstract base class for N-dimensional shapes.

Main implementing classes:

Sphere
box family: box (generator), Box, Cuboid, BaseBox

All geometry objects support batching. Thereby any parameter defining the geometry can be varied along arbitrary batch dims. All batch dimensions are listed in Geometry.shape.

Property getters (@property, such as shape), save for getters, must not depend on any variables marked as variable via __variable_attrs__() as these may be None during tracing. Equality checks must also take this into account.

Expand source code

class Geometry:
    """
    Abstract base class for N-dimensional shapes.

    Main implementing classes:

    * Sphere
    * box family: box (generator), Box, Cuboid, BaseBox

    All geometry objects support batching.
    Thereby any parameter defining the geometry can be varied along arbitrary batch dims.
    All batch dimensions are listed in Geometry.shape.

    Property getters (`@property`, such as `shape`), save for getters, must not depend on any variables marked as *variable* via `__variable_attrs__()` as these may be `None` during tracing.
    Equality checks must also take this into account.
    """

    @property
    def center(self) -> Tensor:
        """
        Center location in single channel dimension.
        """
        raise NotImplementedError(self.__class__)

    @property
    def shape(self) -> Shape:
        """
        The `shape` of a `Geometry` consists of the following dimensions:

        * A single *channel* dimension called `'vector'` specifying the physical space
        * Instance dimensions denote that this geometry consists of multiple copies in the same space
        * Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
        * Batch dimensions indicate non-interacting versions of this geometry for parallelization only.
        """
        raise NotImplementedError(self.__class__)

    @property
    def volume(self) -> Tensor:
        """
        `phi.math.Tensor` representing the volume of each element.
        The result retains batch, spatial and instance dimensions.
        """
        raise NotImplementedError(self.__class__)

    @property
    def faces(self) -> 'Geometry':
        raise NotImplementedError(self.__class__)

    @property
    def face_centers(self) -> Tensor:
        """
        Center of face connecting a pair of cells. Shape `(elements, ~, vector)`.
        Here, `~` represents arbitrary internal dual dimensions, such as `~staggered_direction` or `~elements`.
        Returns 0-vectors for unconnected cells.
        """
        raise NotImplementedError(self.__class__)

    @property
    def face_areas(self) -> Tensor:
        """
        Area of face connecting a pair of cells. Shape `(elements, ~)`.
        Returns 0 for unconnected cells.
        """
        raise NotImplementedError(self.__class__)

    @property
    def face_normals(self) -> Tensor:
        """
        Normal vectors of cell faces, including boundary faces. Shape `(elements, ~, vector)`.
        For meshes, The vectors point out of the primal cells and into the dual cells.

        Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.
        """
        raise NotImplementedError(self.__class__)

    @property
    def boundary_elements(self) -> Dict[str, Dict[str, slice]]:
        """
        Slices on the primal dimensions to mark boundary elements.
        Grids and meshes have no boundary elements and return `{}`.
        Dynamic graphs can define boundary elements for obstacles and walls.

        Returns:
            Map from `name` to slicing `dict`.
        """
        raise NotImplementedError(self.__class__)

    @property
    def boundary_faces(self) -> Dict[str, Dict[str, slice]]:
        """
        Slices on the dual dimensions to mark boundary faces.

        Regular grids use the keys (dim, is_upper) to identify boundaries.
        Unstructured meshes use string identifiers for the boundaries.
        Dynamic graphs return slices along the dual dimensions.

        Returns:
            Map from `name` to slicing `dict`.
        """
        raise NotImplementedError(self.__class__)

    @property
    def face_shape(self) -> Shape:
        """
        Returns:
            Full Shape to identify each face of this `Geometry`, including instance/spatial dimensions for the elements and dual dimensions listing the faces per element.
            If this `Geometry` has no faces, returns an empty `Shape`.
        """
        return None

    @property
    def sets(self) -> Dict[str, Shape]:
        if self.face_shape and self.face_shape != self.shape and self.face_shape.volume > 0:
            return {'center': non_batch(self)-'vector', 'face': self.face_shape.non_batch}
        else:
            return {'center': non_batch(self)-'vector'}

    def get_points(self, set_key: str) -> Tensor:
        if set_key == 'center':
            return self.center
        elif set_key == 'face':
            return self.face_centers
        else:
            raise ValueError(f"Unknown set: '{set_key}'")

    def get_boundary(self, set_key: str) -> Dict[str, Dict[str, slice]]:
        if set_key == 'center':
            return self.boundary_elements
        elif set_key == 'face':
            return self.boundary_faces
        else:
            raise ValueError(f"Unknown set: '{set_key}'")

    @property
    def corners(self) -> Tensor:
        """
        Returns:
            Corner locations as `phiml.math.Tensor`.
            Corners belonging to one object or cell are listed along dual dimensions.
            If the object has no corners, a size-0 tensor with the correct vector and instance dims is returned.
        """
        raise NotImplementedError(self.__class__)

    def integrate_surface(self, face_values: Tensor, divide_volume=False) -> Tensor:
        """
        Multiplies `values´ by the corresponding face area, computes the sum over all faces and divides by the cell volume.
        ∑ values * A.

        Args:
            face_values: Values sampled at the face centers.
            divide_volume: Whether to divide by the cell `volume´

        Returns:
            `Tensor` of values sampled at the centroids.
        """
        result = math.sum(face_values * self.face_areas, self.face_shape.dual)
        return result / self.volume if divide_volume else result

    def integrate_flux(self, flux: Tensor, divide_volume=False) -> Tensor:
        assert 'vector' in flux.shape, f"flux must have a 'vector' dimension but got {flux.shape}"
        result = math.sum(flux.vector @ (self.face_normals * self.face_areas).vector, self.face_shape.dual)
        return result / self.volume if divide_volume else result

    # def resample_to_faces(self, values: Tensor, boundary: Extrapolation, **kwargs):
    #     raise NotImplementedError(self.__class__)
    #
    # def resample_to_centers(self, values: Tensor, boundary: Extrapolation, **kwargs):
    #     raise NotImplementedError(self.__class__)
    #
    # def centered_gradient_of(self, values: Tensor, boundary: Extrapolation, dims=None, **kwargs):
    #     raise NotImplementedError(self.__class__)
    #
    # def staggered_gradient_of(self, values: Tensor, boundary: Extrapolation, dims=None, **kwargs):
    #     raise NotImplementedError(self.__class__)
    #
    # def divergence_of(self, values: Tensor, boundary: Extrapolation, dims=None, **kwargs):
    #     raise NotImplementedError(self.__class__)
    #
    # def laplace_of(self, values: Tensor, boundary: Extrapolation, dims=None, **kwargs):
    #     raise NotImplementedError(self.__class__)
    #
    # def centered_curl_of(self, values: Tensor, boundary: Extrapolation, dims=None, **kwargs):
    #     raise NotImplementedError(self.__class__)
    #
    # def staggered_curl_of(self, values: Tensor, boundary: Extrapolation, dims=None, **kwargs):
    #     raise NotImplementedError(self.__class__)

    def unstack(self, dimension: str) -> tuple:
        """
        Unstacks this Geometry along the given dimension.
        The shapes of the returned geometries are reduced by `dimension`.

        Args:
            dimension: dimension along which to unstack

        Returns:
            geometries: tuple of length equal to `geometry.shape.get_size(dimension)`
        """
        warnings.warn(f"Geometry.unstack() is deprecated. Use math.unstack(geometry) instead.", DeprecationWarning)
        return math.unstack(self, dimension)

    @property
    def spatial_rank(self) -> int:
        """ Number of spatial dimensions of the geometry, 1 = 1D, 2 = 2D, 3 = 3D, etc. """
        return self.shape.get_size('vector')

    def lies_inside(self, location: Tensor) -> Tensor:
        """
        Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

        When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

        Args:
          location: float tensor of shape (batch_size, ..., rank)

        Returns:
          bool tensor of shape (*location.shape[:-1], 1).

        """
        raise NotImplementedError(self.__class__)

    def approximate_closest_surface(self, location: Tensor) -> Tuple[Tensor, Tensor, Tensor, Tensor, Tensor]:
        """
        Find the closest surface face of this geometry given a point that can be outside or inside the geometry.

        Args:
            location: `Tensor` with a single channel dimension called vector. Can have arbitrary other dimensions.

        Returns:
            signed_distance: Scalar signed distance from `location`  to the closest point on the surface.
                Positive values indicate the point lies outside the geometry, negative values indicate the point lies inside the geometry.
            delta: Vector-valued distance vector from `location` to the closest point on the surface.
            normal: Closest surface normal vector.
            offset: Min distance of a surface-tangential plane from 0 as a scalar.
            face_index: (Optional) An index vector pointing at the closest face.
        """
        raise NotImplementedError(self.__class__)

    def approximate_signed_distance(self, location: Tensor) -> Tensor:
        """
        Computes the approximate distance from location to the surface of the geometry.
        Locations outside return positive values, inside negative values and zero exactly at the boundary.

        The exact distance metric used depends on the geometry.
        The approximation holds close to the surface and the distance grows to infinity as the location is moved infinitely far from the geometry.
        The distance metric is differentiable and its gradients are bounded at every point in space.

        When dealing with unions or collections of geometries (instance dimensions), the shortest distance to any instance is returned.
        This also holds for negative distances.

        Args:
            location: `Tensor` with one channel dim `vector` matching the geometry's `vector` dim.

        Returns:
            Float `Tensor`
        """
        raise NotImplementedError(self.__class__)

    def approximate_fraction_inside(self, other_geometry: 'Geometry', balance: Union[Tensor, Number] = 0.5) -> Tensor:
        """
        Computes the approximate overlap between the geometry and a small other geometry.
        Returns 1.0 if `other_geometry` is fully enclosed in this geometry and 0.0 if there is no overlap.
        Close to the surface of this geometry, the fraction filled is differentiable w.r.t. the location and size of `other_geometry`.

        To call this method on batches of geometries of same shape, pass a batched Geometry instance.
        The result tensor will match the batch shape of `other_geometry`.

        The result may only be accurate in special cases.
        The given geometries may be approximated as spheres or boxes using `bounding_radius()` and `bounding_half_extent()`.

        The default implementation of this method approximates other_geometry as a Sphere and computes the fraction using `approximate_signed_distance()`.

        Args:
            other_geometry: `Geometry` or geometry batch for which to compute the overlap with `self`.
            balance: Mid-level between 0 and 1, default 0.5.
                This value is returned when exactly half of `other_geometry` lies inside `self`.
                `0.5 < balance <= 1` makes `self` seem larger while `0 <= balance < 0.5`makes `self` seem smaller.

        Returns:
          fraction of cell volume lying inside the geometry. float tensor of shape (other_geometry.batch_shape, 1).

        """
        assert isinstance(other_geometry, Geometry)
        radius = other_geometry.bounding_radius()
        location = other_geometry.center
        distance = self.approximate_signed_distance(location)
        inside_fraction = balance - distance / radius
        inside_fraction = math.clip(inside_fraction, 0, 1)
        return inside_fraction

    def push(self, positions: Tensor, outward: bool = True, shift_amount: float = 0) -> Tensor:
        """
        Shifts positions either into or out of geometry.

        Args:
            positions: Tensor holding positions to shift
            outward: Flag for indicating inward (False) or outward (True) shift
            shift_amount: Minimum distance between positions and surface after shifting.

        Returns:
            Tensor holding shifted positions.
        """
        from ._geom_ops import expel
        return expel(self, positions, min_separation=shift_amount, invert=not outward)

    def sample_uniform(self, *shape: math.Shape) -> Tensor:
        """
        Samples uniformly distributed random points inside this volume.

        Args:
            *shape: How many points to sample per individual geometry.

        Returns:
            `Tensor` containing all dimensions from `Geometry.shape`, `shape` as well as a `channel` dimension `vector` matching the dimensionality of this `Geometry`.
        """
        raise NotImplementedError(self.__class__)

    def bounding_radius(self) -> Tensor:
        """
        Returns the radius of a Sphere object that fully encloses this geometry.
        The sphere is centered at the center of this geometry.

        If this geometry consists of multiple parts listed along instance/spatial dims, these dims are retained, giving the bounds of each part.
        If these dims are not present on the result, all parts are assumed to have the same bounds.
        """
        raise NotImplementedError(self.__class__)

    def bounding_half_extent(self) -> Tensor:
        """
        The bounding half-extent sets a limit on the outer-most point for each coordinate axis.
        Each component is non-negative.

        Let the bounding half-extent have value `e` in dimension `d` (`extent[...,d] = e`).
        Then, no point of the geometry lies further away from its center point than `e` along `d` (in both axis directions).

        If this geometry consists of multiple parts listed along instance/spatial dims, these dims are retained, giving the bounds of each part.
        If these dims are not present on the result, all parts are assumed to have the same bounds.
        """
        raise NotImplementedError(self.__class__)

    def bounding_box(self) -> 'BaseBox':
        """
        Returns the approximately smallest axis-aligned box that contains this `Geometry`.
        The center of the box may not be equal to `self.center`.

        Returns:
            `Box` or `Cuboid` that fully contains this `Geometry`.
        """
        center = self.center
        half = self.bounding_half_extent()
        min_vec = math.min(center - half, dim=center.shape.non_batch.non_channel)
        max_vec = math.max(center + half, dim=center.shape.non_batch.non_channel)
        from ._box import Box
        return Box(min_vec, max_vec)

    def shifted(self, delta: Tensor) -> 'Geometry':
        """
        Returns a translated version of this geometry.

        See Also:
            `Geometry.at()`.

        Args:
          delta: direction vector
          delta: Tensor:

        Returns:
          Geometry: shifted geometry

        """
        return self.at(self.center + delta)

    def at(self, center: Tensor) -> 'Geometry':
        """
        Returns a copy of this `Geometry` with the center at `center`.
        This is equal to calling `self @ center`.

        See Also:
            `Geometry.shifted()`.

        Args:
            center: New center as `Tensor`.

        Returns:
            `Geometry`.
        """
        raise NotImplementedError(self.__class__)

    def __matmul__(self, other):
        if isinstance(other, (Tensor, float, int)):
            return self.at(other)
        return NotImplemented

    def rotated(self, angle: Union[float, Tensor]) -> 'Geometry':
        """
        Returns a rotated version of this geometry.
        The geometry is rotated about its center point.

        Args:
            angle: Delta rotation.
                Either

                * Angle(s): scalar angle in 2d or euler angles along `vector` in 3D or higher.
                * Matrix: d⨯d rotation matrix

        Returns:
            Rotated `Geometry`
        """
        raise NotImplementedError(self.__class__)

    def scaled(self, factor: Union[float, Tensor]) -> 'Geometry':
        """
        Scales each individual geometry by `factor`.
        The individual `center` points act as pivots for the operation.

        Args:
            factor:

        Returns:

        """
        raise NotImplementedError(self.__class__)

    def __invert__(self):
        return InvertedGeometry(self)

    def __eq__(self, other):
        """
        Slow equality check.
        Unlike `==`, this method compares all tensor elements to check whether they are equal.
        Use `==` for a faster check which only checks whether the referenced tensors are the same.

        See Also:
            `shallow_equals()`
        """
        def tensor_equality(a, b):
            if a is None or b is None:
                return True  # stored mode, tensors unavailable
            return math.close(a, b, rel_tolerance=1e-5, equal_nan=True)
        differences = find_differences(self, other, attr_type=variable_attributes, tensor_equality=tensor_equality)
        return not differences

    def shallow_equals(self, other):
        """
        Quick equality check.
        May return `False` even if `other == self`.
        However, if `True` is returned, the geometries are guaranteed to be equal.

        The `shallow_equals()` check does not compare all tensor elements but merely checks whether the same tensors are referenced.
        """
        differences = find_differences(self, other, compare_tensors_by_id=True)
        return not differences

    @staticmethod
    def __stack__(values: tuple, dim: Shape, **kwargs) -> 'Geometry':
        if all(type(v) == type(values[0]) for v in values):
            return NotImplemented  # let attributes be stacked
        else:
            from ._geom_ops import GeometryStack
            set_op = kwargs.get('set_op')
            return GeometryStack(math.layout(values, dim), set_op)

    def __flatten__(self, flat_dim: Shape, flatten_batch: bool, **kwargs) -> 'Geometry':
        dims = self.shape.without('vector')
        if not flatten_batch:
            dims = dims.non_batch
        return math.pack_dims(self, dims, flat_dim, **kwargs)

    def __ne__(self, other):
        return not self == other

    def __hash__(self):
        return id(self.__class__) + hash(self.shape)

    def __repr__(self):
        return f"{self.__class__.__name__}{self.shape}"

    def __getitem__(self, item):
        raise NotImplementedError
        # attrs = {a: getattr(self, a)[item] for a in variable_attributes(self)}
        # return copy_with(self, **attrs)

    def __getattr__(self, name: str) -> BoundDim:
        return BoundDim(self, name)

Subclasses

phi.geom._box.BaseBox
phi.geom._cylinder.Cylinder
phi.geom._embed._EmbeddedGeometry
phi.geom._geom.InvertedGeometry
phi.geom._geom.NoGeometry
phi.geom._geom.Point
phi.geom._geom_ops.GeometryStack
phi.geom._graph.Graph
phi.geom._heightmap.Heightmap
phi.geom._mesh.Mesh
phi.geom._sdf.SDF
phi.geom._sdf_grid.SDFGrid
phi.geom._sphere.Sphere

Instance variables

prop boundary_elements : Dict[str, Dict[str, slice]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[str, Dict[str, slice]]:
    """
    Slices on the primal dimensions to mark boundary elements.
    Grids and meshes have no boundary elements and return `{}`.
    Dynamic graphs can define boundary elements for obstacles and walls.

    Returns:
        Map from `name` to slicing `dict`.
    """
    raise NotImplementedError(self.__class__)

prop boundary_faces : Dict[str, Dict[str, slice]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[str, Dict[str, slice]]:
    """
    Slices on the dual dimensions to mark boundary faces.

    Regular grids use the keys (dim, is_upper) to identify boundaries.
    Unstructured meshes use string identifiers for the boundaries.
    Dynamic graphs return slices along the dual dimensions.

    Returns:
        Map from `name` to slicing `dict`.
    """
    raise NotImplementedError(self.__class__)

prop center : phiml.math._tensors.Tensor

Center location in single channel dimension.

Expand source code

@property
def center(self) -> Tensor:
    """
    Center location in single channel dimension.
    """
    raise NotImplementedError(self.__class__)

prop corners : phiml.math._tensors.Tensor

Returns

Expand source code

@property
def corners(self) -> Tensor:
    """
    Returns:
        Corner locations as `phiml.math.Tensor`.
        Corners belonging to one object or cell are listed along dual dimensions.
        If the object has no corners, a size-0 tensor with the correct vector and instance dims is returned.
    """
    raise NotImplementedError(self.__class__)

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    """
    Area of face connecting a pair of cells. Shape `(elements, ~)`.
    Returns 0 for unconnected cells.
    """
    raise NotImplementedError(self.__class__)

prop face_centers : phiml.math._tensors.Tensor

Expand source code

@property
def face_centers(self) -> Tensor:
    """
    Center of face connecting a pair of cells. Shape `(elements, ~, vector)`.
    Here, `~` represents arbitrary internal dual dimensions, such as `~staggered_direction` or `~elements`.
    Returns 0-vectors for unconnected cells.
    """
    raise NotImplementedError(self.__class__)

prop face_normals : phiml.math._tensors.Tensor

Normal vectors of cell faces, including boundary faces. Shape (elements, ~, vector). For meshes, The vectors point out of the primal cells and into the dual cells.

Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.

Expand source code

@property
def face_normals(self) -> Tensor:
    """
    Normal vectors of cell faces, including boundary faces. Shape `(elements, ~, vector)`.
    For meshes, The vectors point out of the primal cells and into the dual cells.

    Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.
    """
    raise NotImplementedError(self.__class__)

prop face_shape : phiml.math._shape.Shape

Returns

Expand source code

@property
def face_shape(self) -> Shape:
    """
    Returns:
        Full Shape to identify each face of this `Geometry`, including instance/spatial dimensions for the elements and dual dimensions listing the faces per element.
        If this `Geometry` has no faces, returns an empty `Shape`.
    """
    return None

prop faces : Geometry

Expand source code

@property
def faces(self) -> 'Geometry':
    raise NotImplementedError(self.__class__)

prop sets : Dict[str, phiml.math._shape.Shape]

Expand source code

@property
def sets(self) -> Dict[str, Shape]:
    if self.face_shape and self.face_shape != self.shape and self.face_shape.volume > 0:
        return {'center': non_batch(self)-'vector', 'face': self.face_shape.non_batch}
    else:
        return {'center': non_batch(self)-'vector'}

prop shape : phiml.math._shape.Shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self) -> Shape:
    """
    The `shape` of a `Geometry` consists of the following dimensions:

    * A single *channel* dimension called `'vector'` specifying the physical space
    * Instance dimensions denote that this geometry consists of multiple copies in the same space
    * Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
    * Batch dimensions indicate non-interacting versions of this geometry for parallelization only.
    """
    raise NotImplementedError(self.__class__)

prop spatial_rank : int

Number of spatial dimensions of the geometry, 1 = 1D, 2 = 2D, 3 = 3D, etc.

Expand source code

@property
def spatial_rank(self) -> int:
    """ Number of spatial dimensions of the geometry, 1 = 1D, 2 = 2D, 3 = 3D, etc. """
    return self.shape.get_size('vector')

prop volume : phiml.math._tensors.Tensor

phi.math.Tensor representing the volume of each element. The result retains batch, spatial and instance dimensions.

Expand source code

@property
def volume(self) -> Tensor:
    """
    `phi.math.Tensor` representing the volume of each element.
    The result retains batch, spatial and instance dimensions.
    """
    raise NotImplementedError(self.__class__)

Methods

Find the closest surface face of this geometry given a point that can be outside or inside the geometry.

Args

location: Tensor with a single channel dimension called vector. Can have arbitrary other dimensions.

Returns

signed_distance: Scalar signed distance from location to the closest point on the surface. Positive values indicate the point lies outside the geometry, negative values indicate the point lies inside the geometry.
delta: Vector-valued distance vector from location to the closest point on the surface.
normal: Closest surface normal vector.
offset: Min distance of a surface-tangential plane from 0 as a scalar.
face_index: (Optional) An index vector pointing at the closest face.

def approximate_fraction_inside(self, other_geometry: Geometry, balance: phiml.math._tensors.Tensor | numbers.Number = 0.5) ‑> phiml.math._tensors.Tensor

Computes the approximate overlap between the geometry and a small other geometry. Returns 1.0 if other_geometry is fully enclosed in this geometry and 0.0 if there is no overlap. Close to the surface of this geometry, the fraction filled is differentiable w.r.t. the location and size of other_geometry.

To call this method on batches of geometries of same shape, pass a batched Geometry instance. The result tensor will match the batch shape of other_geometry.

The result may only be accurate in special cases. The given geometries may be approximated as spheres or boxes using bounding_radius() and bounding_half_extent().

The default implementation of this method approximates other_geometry as a Sphere and computes the fraction using approximate_signed_distance().

Args

other_geometry: Geometry or geometry batch for which to compute the overlap with self.
balance: Mid-level between 0 and 1, default 0.5. This value is returned when exactly half of other_geometry lies inside self. 0.5 < balance <= 1 makes self seem larger while 0 <= balance < 0.5makes self seem smaller.

Returns

fraction of cell volume lying inside the geometry. float tensor of shape (other_geometry.batch_shape, 1).

def approximate_signed_distance(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Computes the approximate distance from location to the surface of the geometry. Locations outside return positive values, inside negative values and zero exactly at the boundary.

When dealing with unions or collections of geometries (instance dimensions), the shortest distance to any instance is returned. This also holds for negative distances.

Args

location: Tensor with one channel dim vector matching the geometry's vector dim.

Returns

Float Tensor

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_box(self)

Returns the approximately smallest axis-aligned box that contains this Geometry. The center of the box may not be equal to self.center.

Returns

Box or Cuboid that fully contains this Geometry.

def bounding_half_extent(self) ‑> phiml.math._tensors.Tensor

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def bounding_radius(self) ‑> phiml.math._tensors.Tensor

Returns the radius of a Sphere object that fully encloses this geometry. The sphere is centered at the center of this geometry.

def get_boundary(self, set_key: str) ‑> Dict[str, Dict[str, slice]]

def get_points(self, set_key: str) ‑> phiml.math._tensors.Tensor

def integrate_flux(self, flux: phiml.math._tensors.Tensor, divide_volume=False) ‑> phiml.math._tensors.Tensor

def integrate_surface(self, face_values: phiml.math._tensors.Tensor, divide_volume=False) ‑> phiml.math._tensors.Tensor

Multiplies `values´ by the corresponding face area, computes the sum over all faces and divides by the cell volume. ∑ values * A.

Args

face_values: Values sampled at the face centers.
divide_volume: Whether to divide by the cell `volume´

Returns

Tensor of values sampled at the centroids.

def lies_inside(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def push(self, positions: phiml.math._tensors.Tensor, outward: bool = True, shift_amount: float = 0) ‑> phiml.math._tensors.Tensor

Shifts positions either into or out of geometry.

Args

positions: Tensor holding positions to shift
outward: Flag for indicating inward (False) or outward (True) shift
shift_amount: Minimum distance between positions and surface after shifting.

Returns

Tensor holding shifted positions.

def rotated(self, angle: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def sample_uniform(self, *shape: phiml.math._shape.Shape) ‑> phiml.math._tensors.Tensor

Samples uniformly distributed random points inside this volume.

Args

*shape: How many points to sample per individual geometry.

Returns

Tensor containing all dimensions from Geometry.shape, shape as well as a channel dimension vector matching the dimensionality of this Geometry.

def scaled(self, factor: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Scales each individual geometry by factor. The individual center points act as pivots for the operation.

Args

factor: Returns:

def shallow_equals(self, other)

Quick equality check. May return False even if other == self. However, if True is returned, the geometries are guaranteed to be equal.

The shallow_equals() check does not compare all tensor elements but merely checks whether the same tensors are referenced.

def shifted(self, delta: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a translated version of this geometry.

See Also: Geometry.at().

Args

delta: direction vector
delta: Tensor:

Returns

Geometry: shifted geometry

def unstack(self, dimension: str) ‑> tuple

Unstacks this Geometry along the given dimension. The shapes of the returned geometries are reduced by dimension.

Args

dimension: dimension along which to unstack

Returns

geometries: tuple of length equal to geometry.shape.get_size(dimension)

class GeometryException (*args, **kwargs)

Raised when an operation is fundamentally not possible for a Geometry. Possible causes:

Trying to get the interior of a non-surface Geometry
Trying to get the surface of a point-like Geometry

Expand source code

class GeometryException(BaseException):
    """
    Raised when an operation is fundamentally not possible for a `Geometry`.
    Possible causes:

    * Trying to get the interior of a non-surface `Geometry`
    * Trying to get the surface of a point-like `Geometry`
    """

Ancestors

builtins.BaseException

class Graph (nodes: phi.geom._geom.Geometry | phiml.math._tensors.Tensor, edges: phiml.math._tensors.Tensor, boundary: Dict[str, Dict[str, slice]], deltas: phiml.math._tensors.Tensor | None = None, distances: phiml.math._tensors.Tensor | None = None, bounding_distance: float | phiml.math._tensors.Tensor | None = None)

A graph consists of multiple geometry nodes and corresponding edge information.

Edges are stored as a Tensor with the same axes ad geometry plus their dual counterparts. Additional dimensions can be added to edges to store vector-valued connectivity weights.

Create a graph where nodes are connected by edges.

Args

nodes: Geometry collection or Tensor to denote points.
edges: Edge weight matrix. Must have the instance and spatial dims of nodes plus their dual counterparts.
boundary: Marks ranges of nodes as boundary elements.
deltas: (Optional) Pre-computed position difference matrix.
distances: (Optional) Pre-computed distance matrix.
bounding_distance: (Optional) Pre-computed distance bounds. No distance is larger than this value. If True, will be computed now, if False, will not be computed.

Expand source code

class Graph(Geometry):
    """
    A graph consists of multiple geometry nodes and corresponding edge information.

    Edges are stored as a Tensor with the same axes ad `geometry` plus their dual counterparts.
    Additional dimensions can be added to `edges` to store vector-valued connectivity weights.
    """

    def __init__(self,
                 nodes: Union[Geometry, Tensor],
                 edges: Tensor,
                 boundary: Dict[str, Dict[str, slice]],
                 deltas: Optional[Tensor] = None,
                 distances: Optional[Tensor] = None,
                 bounding_distance: Tensor | float | None = None):
        """
        Create a graph where `nodes` are connected by `edges`.

        Args:
            nodes: `Geometry` collection or `Tensor` to denote points.
            edges: Edge weight matrix. Must have the instance and spatial dims of `nodes` plus their dual counterparts.
            boundary: Marks ranges of nodes as boundary elements.
            deltas: (Optional) Pre-computed position difference matrix.
            distances: (Optional) Pre-computed distance matrix.
            bounding_distance: (Optional) Pre-computed distance bounds. No distance is larger than this value. If `True`, will be computed now, if `False`, will not be computed.
        """
        assert isinstance(nodes, Geometry), f"nodes must be a Geometry  but got {nodes}"
        node_dims = non_batch(nodes).non_channel
        assert node_dims in edges.shape and edges.shape.dual.rank == node_dims.rank, f"edges must contain all node dims {node_dims} as primal and dual but got {edges.shape}"
        self._nodes: Geometry = nodes if isinstance(nodes, Geometry) else Point(nodes)
        self._edges = edges
        self._boundary = boundary
        self._deltas = deltas
        self._distances = distances
        self._connectivity = math.tensor_like(edges, 1) if math.is_sparse(edges) else (edges != 0) & ~math.is_nan(edges)
        if isinstance(bounding_distance, bool):
            self._bounding_distance = math.max(self._distances) if bounding_distance else None
        else:
            self._bounding_distance = bounding_distance

    def __variable_attrs__(self):
        return '_nodes', '_edges', '_deltas', '_distances', '_connectivity'

    def __value_attrs__(self):
        return '_nodes',

    @property
    def edges(self):
        return self._edges

    @property
    def connectivity(self) -> Tensor:
        return self._connectivity

    @property
    def nodes(self) -> Geometry:
        return self._nodes

    def as_points(self):
        return Graph(Point(self._nodes.center), self._edges, self._boundary, self._deltas, self._distances, self._bounding_distance)

    @property
    def deltas(self):
        return self._deltas

    @property
    def unit_deltas(self):
        return math.safe_div(self._deltas, self._distances)

    @property
    def distances(self):
        return self._distances

    @property
    def bounding_distance(self) -> Optional[Tensor]:
        return self._bounding_distance

    @property
    def center(self) -> Tensor:
        return self._nodes.center

    @property
    def shape(self) -> Shape:
        return self._nodes.shape

    @property
    def volume(self) -> Tensor:
        return self._nodes.volume

    @property
    def faces(self) -> 'Geometry':
        raise NotImplementedError

    @property
    def face_centers(self) -> Tensor:
        raise NotImplementedError

    @property
    def face_areas(self) -> Tensor:
        raise NotImplementedError

    @property
    def face_normals(self) -> Tensor:
        raise NotImplementedError

    @property
    def boundary_elements(self) -> Dict[str, Dict[str, slice]]:
        return self._boundary

    @property
    def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
        raise NotImplementedError  # connections between boundary elements

    @property
    def face_shape(self) -> Shape:
        return non_batch(self._edges).non_channel

    def lies_inside(self, location: Tensor) -> Tensor:
        raise NotImplementedError

    def approximate_closest_surface(self, location: Tensor) -> Tuple[Tensor, Tensor, Tensor, Tensor, Tensor]:
        raise NotImplementedError

    def approximate_signed_distance(self, location: Tensor) -> Tensor:
        raise NotImplementedError

    def sample_uniform(self, *shape: math.Shape) -> Tensor:
        raise NotImplementedError

    def bounding_radius(self) -> Tensor:
        return self._nodes.bounding_radius()

    def bounding_half_extent(self) -> Tensor:
        return self._nodes.bounding_half_extent()

    def at(self, center: Tensor) -> 'Geometry':
        raise NotImplementedError("Changing the node positions of a Graph is not supported as it would invalidate distances.")
        # warnings.warn("Changing the node positions of a graph triggers re-evaluation of distances.", RuntimeWarning, stacklevel=2)
        # return Graph(self.nodes.at(center), self._edges, self._boundary, bounding_distance=self._bounding_distance is not None)

    def shifted(self, delta: Tensor) -> 'Geometry':
        if non_batch(delta).non_channel.only(self._nodes.shape) and self._deltas is not None:  # shift varies between elements
            raise NotImplementedError("Shifting the node positions of a Graph is not supported as it would invalidate distances.")
        return Graph(self.nodes.shifted(delta), self._edges, self._boundary, deltas=self._deltas, distances=self._distances, bounding_distance=self._bounding_distance is not None)

    def rotated(self, angle: Union[float, Tensor]) -> 'Geometry':
        raise NotImplementedError

    def scaled(self, factor: Union[float, Tensor]) -> 'Geometry':
        raise NotImplementedError

    def __getitem__(self, item):
        item = slicing_dict(self, item)
        node_dims = non_batch(self._nodes).non_channel
        edge_sel = {}
        for i, (dim, sel) in enumerate(item.items()):
            if dim in node_dims:
                dual_dim = '~' + dim
                if dual_dim not in self._edges.shape:
                    dual_dim = dual(self._edges).shape.names[i]
                edge_sel[dim] = edge_sel[dual_dim] = sel
            elif dim in batch(self):
                edge_sel[dim] = sel
        deltas = self._deltas[edge_sel] if self._deltas is not None else None
        distances = self._distances[edge_sel] if self._distances is not None else None
        bounding_distance = self._bounding_distance[item] if self._bounding_distance is not None else None
        return Graph(self._nodes[item], self._edges[edge_sel], self._boundary, deltas, distances, bounding_distance)

Ancestors

phi.geom._geom.Geometry

Instance variables

prop boundary_elements : Dict[str, Dict[str, slice]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[str, Dict[str, slice]]:
    return self._boundary

prop boundary_faces : Dict[Any, Dict[str, slice]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
    raise NotImplementedError  # connections between boundary elements

prop bounding_distance : phiml.math._tensors.Tensor | None

Expand source code

@property
def bounding_distance(self) -> Optional[Tensor]:
    return self._bounding_distance

prop center : phiml.math._tensors.Tensor

Center location in single channel dimension.

Expand source code

@property
def center(self) -> Tensor:
    return self._nodes.center

prop connectivity : phiml.math._tensors.Tensor

Expand source code

@property
def connectivity(self) -> Tensor:
    return self._connectivity

prop deltas

Expand source code

@property
def deltas(self):
    return self._deltas

prop distances

Expand source code

@property
def distances(self):
    return self._distances

prop edges

Expand source code

@property
def edges(self):
    return self._edges

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    raise NotImplementedError

prop face_centers : phiml.math._tensors.Tensor

Expand source code

@property
def face_centers(self) -> Tensor:
    raise NotImplementedError

prop face_normals : phiml.math._tensors.Tensor

Normal vectors of cell faces, including boundary faces. Shape (elements, ~, vector). For meshes, The vectors point out of the primal cells and into the dual cells.

Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.

Expand source code

@property
def face_normals(self) -> Tensor:
    raise NotImplementedError

prop face_shape : phiml.math._shape.Shape

Returns

Expand source code

@property
def face_shape(self) -> Shape:
    return non_batch(self._edges).non_channel

prop faces : Geometry

Expand source code

@property
def faces(self) -> 'Geometry':
    raise NotImplementedError

prop nodes : phi.geom._geom.Geometry

Expand source code

@property
def nodes(self) -> Geometry:
    return self._nodes

prop shape : phiml.math._shape.Shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self) -> Shape:
    return self._nodes.shape

prop unit_deltas

Expand source code

@property
def unit_deltas(self):
    return math.safe_div(self._deltas, self._distances)

prop volume : phiml.math._tensors.Tensor

phi.math.Tensor representing the volume of each element. The result retains batch, spatial and instance dimensions.

Expand source code

@property
def volume(self) -> Tensor:
    return self._nodes.volume

Methods

Find the closest surface face of this geometry given a point that can be outside or inside the geometry.

Args

location: Tensor with a single channel dimension called vector. Can have arbitrary other dimensions.

Returns

signed_distance: Scalar signed distance from location to the closest point on the surface. Positive values indicate the point lies outside the geometry, negative values indicate the point lies inside the geometry.
delta: Vector-valued distance vector from location to the closest point on the surface.
normal: Closest surface normal vector.
offset: Min distance of a surface-tangential plane from 0 as a scalar.
face_index: (Optional) An index vector pointing at the closest face.

def approximate_signed_distance(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Computes the approximate distance from location to the surface of the geometry. Locations outside return positive values, inside negative values and zero exactly at the boundary.

When dealing with unions or collections of geometries (instance dimensions), the shortest distance to any instance is returned. This also holds for negative distances.

Args

location: Tensor with one channel dim vector matching the geometry's vector dim.

Returns

Float Tensor

def as_points(self)

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_half_extent(self) ‑> phiml.math._tensors.Tensor

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def bounding_radius(self) ‑> phiml.math._tensors.Tensor

Returns the radius of a Sphere object that fully encloses this geometry. The sphere is centered at the center of this geometry.

def lies_inside(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def rotated(self, angle: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def sample_uniform(self, *shape: phiml.math._shape.Shape) ‑> phiml.math._tensors.Tensor

Samples uniformly distributed random points inside this volume.

Args

*shape: How many points to sample per individual geometry.

Returns

Tensor containing all dimensions from Geometry.shape, shape as well as a channel dimension vector matching the dimensionality of this Geometry.

def scaled(self, factor: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Scales each individual geometry by factor. The individual center points act as pivots for the operation.

Args

factor: Returns:

def shifted(self, delta: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a translated version of this geometry.

See Also: Geometry.at().

Args

delta: direction vector
delta: Tensor:

Returns

Geometry: shifted geometry

class Heightmap (height: phiml.math._tensors.Tensor, bounds: phi.geom._box.Box, max_dist: phiml.math._tensors.Tensor | float, fill_below: bool | phiml.math._tensors.Tensor = True, extrapolation: float | str | phiml.math.extrapolation.Extrapolation = None, faces=None)

Abstract base class for N-dimensional shapes.

Main implementing classes:

Sphere
box family: box (generator), Box, Cuboid, BaseBox

All geometry objects support batching. Thereby any parameter defining the geometry can be varied along arbitrary batch dims. All batch dimensions are listed in Geometry.shape.

Args

height: Heightmap Tensor of absolute (world-space) height values. Scalar height values on a d-1 dimensional grid.
bounds: d-dimensional bounds. Locations outside bounds' can never lie inside this geometry ifextrapolation is None`. Otherwise, only the height dimension is checked. The grid dimensions of bounds must be finite but the height dimension may be infinite to count all values above/below height as inside.
max_dist: Maximum distance up to which the distance approximations should be valid. This does not affect the number of computations performed to compute the distance. Low values increase accuracy close to the surface but trade off possibly very wrong distances further away.
fill_below: Whether the inside is below or above the height values.
extrapolation: Surface height outside bounds´. Can be any validphiml.math.Extrapolation`, such as a constant. If not None, values outside bounds will be checked against the extrapolated height values. Otherwise, values outside bounds always lie on the outside.

Expand source code

class Heightmap(Geometry):

    def __init__(self,
                 height: Tensor,
                 bounds: Box,
                 max_dist: Union[float, Tensor],
                 fill_below: Union[bool, Tensor] = True,
                 extrapolation: Union[float, str, math.Extrapolation] = None,
                 faces=None):
        """

        Args:
            height: Heightmap `Tensor` of absolute (world-space) height values.
                Scalar height values on a d-1 dimensional grid.
            bounds: d-dimensional bounds.
                Locations outside `bounds' can never lie inside this geometry if `extrapolation is None`.
                Otherwise, only the height dimension is checked.
                The grid dimensions of `bounds` must be finite but the height dimension may be infinite to count all values above/below `height` as inside.
            max_dist: Maximum distance up to which the distance approximations should be valid.
                This does not affect the number of computations performed to compute the distance.
                Low values increase accuracy close to the surface but trade off possibly very wrong distances further away.
            fill_below: Whether the inside is below or above the height values.
            extrapolation: Surface height outside `bounds´. Can be any valid `phiml.math.Extrapolation`, such as a constant.
                If not `None`, values outside `bounds` will be checked against the extrapolated `height` values.
                Otherwise, values outside `bounds` always lie on the outside.
        """
        assert channel(height).is_empty, f"height must be a scalar quantity but got {height.shape}"
        assert spatial(height), f"height field must have at least one spatial dim but got {height}"
        assert bounds.vector.size == spatial(height).rank + 1, f"bounds must include the spatial grid dimensions {spatial(height)} and the height dimension but got {bounds}"
        dims = bounds.vector.item_names
        self._hdim = spatial(*dims).without(height.shape).name
        if math.all_available(height, bounds.lower, bounds.upper):
            assert bounds[self._hdim].lies_inside(height).all, f"All height values should be within the {self._hdim}-range given by bounds but got height={height}"
        self._height = height
        self._fill_below = wrap(fill_below)
        self._bounds = bounds
        self._max_dist = wrap(max_dist)
        self._extrapolation = math.as_extrapolation(extrapolation)
        if faces is None:
            proj_faces = build_faces(self)
            with numpy.errstate(divide='ignore', invalid='ignore'):
                secondary_idx = math.map(find_most_important_neighbor, proj_faces, self.dx, self.resolution, self._hdim, self._fill_below, self._max_dist, dims=instance, unwrap_scalars=False)
                secondary_faces = math.map(math.gather, proj_faces, secondary_idx, dims=instance)
            self._faces: Face = stack([proj_faces, *unstack(secondary_faces, 'side')], batch(consider='self,outside,inside'), expand_values=True)
            self._faces = cached(self._faces)  # otherwise, this may get expanded during tracing
        else:
            self._faces = faces

    @property
    def height(self):
        return self._height

    @property
    def bounds(self):
        return self._bounds

    @property
    def max_dist(self):
        return self._max_dist

    @property
    def fill_below(self):
        return self._fill_below

    @property
    def extrapolation(self):
        return self._extrapolation

    @property
    def shape(self) -> Shape:
        return (self._height.shape - 1) & channel(self._bounds)

    @property
    def resolution(self):
        return spatial(self._height) - 1

    @property
    def grid_bounds(self):
        return self._bounds[self.resolution.name_list]

    @property
    def up(self):
        dims = self._bounds.vector.item_names
        height_unit = vec(**{d: 1 if d == self._hdim else 0 for d in dims})
        return math.where(self._fill_below, height_unit, -height_unit)

    @property
    def dx(self):
        return self._bounds.size[self.resolution.name_list] / spatial(self.resolution)

    @property
    def vertices(self):
        hdim = self._hdim
        space = self.vector.item_names
        pos = self.grid_bounds.local_to_global(math.meshgrid(spatial(self._height)) / self.resolution)
        vert = stack({dim: self.height if dim == hdim else pos[dim] for dim in space}, channel('vector'))
        return vert

    def lies_inside(self, location: Tensor) -> Tensor:
        location = rename_dims(location, self.resolution.names, ['loc_' + n for n in self.resolution.names])
        projected_loc = location[self.resolution.name_list]
        @math.map_i2b
        def lies_inside_(height, grid_bounds, bounds, fill_below, extrapolation):
            float_idx = (projected_loc - grid_bounds.lower) / grid_bounds.size * self.resolution
            if extrapolation is None:
                within_bounds = bounds.lies_inside(location)
            else:
                within_bounds = bounds[self._hdim].lies_inside(location[self._hdim])
            surface_height = math.grid_sample(height, float_idx - 1, math.NAN if extrapolation is None else extrapolation)
            is_below = location[self._hdim] <= surface_height
            inside = is_below == fill_below
            result = math.where(within_bounds, inside, False)
            return rename_dims(result, ['loc_' + n for n in self.resolution.names], self.resolution.names)
        return math.any(lies_inside_(self._height, self.grid_bounds, self._bounds, self._fill_below, self._extrapolation), instance(self))

    def approximate_closest_surface(self, location: Tensor) -> Tuple[Tensor, Tensor, Tensor, Tensor, Tensor]:
        grid_bounds = math.i2b(self.grid_bounds)
        faces = math.i2b(self._faces)
        cell_idx = cell_index(location, grid_bounds, self.resolution, clip=True)
        # --- gather face infos at projected cell ---
        normals = faces.normal[cell_idx]
        offsets = faces.origin_distance[cell_idx]
        face_idx = faces.index[cell_idx]
        # --- test location against all considered faces and boundaries ---
        # distances = plane_sgn_dist(-offsets, normals, location)  # offset has the - convention here
        distances = normals.vector @ location.vector + offsets
        projected_onto_face = location - normals * distances
        projected_idx = cell_index(projected_onto_face, grid_bounds, self.resolution, clip=False)
        projects_onto_face = math.all(projected_idx == face_idx, channel)
        proj_delta = normals * -distances
        # --- if not projected onto face, use distance to highest point instead ---
        delta_highest = faces.extrema_points[cell_idx] - location
        flat_normal = math.vec_normalize(normals[self.resolution.name_list], epsilon=1e-5)
        delta_edge = flat_normal * (delta_highest[self.resolution].vector @ flat_normal.vector)  # project onto flat normal
        delta_edge = concat([delta_edge, delta_highest[[self._hdim]]], 'vector')
        distance_edge = math.vec_length(delta_edge, eps=1e-5)
        delta_highest, distance_edge = math.at_min((delta_highest, distance_edge), distance_edge, 'extremum')
        distance_edge = math.where(distances < 0, -distance_edge, distance_edge)  # copy sign of distances onto distance_edges to always return the signed distance
        distances = math.where(projects_onto_face, distances, distance_edge)
        # --- use closest face from considered ---
        delta = math.where(projects_onto_face, proj_delta, delta_highest)
        return math.at_min((distances, delta, normals, offsets, face_idx), key=abs(distances), dim=batch('consider') & instance(self).as_batch())

    def shallow_equals(self, other):
        return self == other

    def __repr__(self):
        return f"Heightmap {self.resolution}, bounds={self._bounds}"

    def __variable_attrs__(self):
        return '_height', '_bounds', '_max_dist', '_fill_below', '_extrapolation', '_faces'

    def __value_attrs__(self):
        return ()

    def __getitem__(self, item):
        item = slicing_dict(self, item)
        return Heightmap(self._height[item], self._bounds[item], self._max_dist[item], self._fill_below[item], self._extrapolation[item] if self._extrapolation is not None else None, math.slice(self._faces, item))

    def bounding_half_extent(self) -> Tensor:
        h_min, h_max = self._faces.extrema_points[{'consider': 0, 'vector': self._hdim}].extremum
        dh = h_max - h_min
        return stack({d: self.dx[d] if d in self.resolution else dh for d in self.vector.item_names}, channel('vector'), expand_values=True) * .5

    @property
    def center(self) -> Tensor:
        return self._faces.center.consider[0]

    @property
    def volume(self) -> Tensor:
        return math.prod(self.bounding_half_extent() * 2, channel)

    @property
    def faces(self) -> 'Geometry':
        raise NotImplementedError

    @property
    def face_centers(self) -> Tensor:
        return self._faces.center

    @property
    def face_areas(self) -> Tensor:
        raise NotImplementedError

    @property
    def face_normals(self) -> Tensor:
        return self._faces.normal

    @property
    def boundary_elements(self) -> Dict[Any, Dict[str, slice]]:
        return {}

    @property
    def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
        return {}

    @property
    def face_shape(self) -> Shape:
        return non_channel(self._faces.center)

    def approximate_signed_distance(self, location: Tensor) -> Tensor:
        return self.approximate_closest_surface(location)[0]

    def sample_uniform(self, *shape: math.Shape) -> Tensor:
        raise NotImplementedError

    def bounding_radius(self) -> Tensor:
        return self._bounds.bounding_radius()

    def at(self, center: Tensor) -> 'Geometry':
        raise NotImplementedError

    def rotated(self, angle: Union[float, Tensor]) -> 'Geometry':
        raise NotImplementedError

    def scaled(self, factor: Union[float, Tensor]) -> 'Geometry':
        raise NotImplementedError

Ancestors

phi.geom._geom.Geometry

Instance variables

prop boundary_elements : Dict[Any, Dict[str, slice]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[Any, Dict[str, slice]]:
    return {}

prop boundary_faces : Dict[Any, Dict[str, slice]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
    return {}

prop bounds

Expand source code

@property
def bounds(self):
    return self._bounds

prop center : phiml.math._tensors.Tensor

Center location in single channel dimension.

Expand source code

@property
def center(self) -> Tensor:
    return self._faces.center.consider[0]

prop dx

Expand source code

@property
def dx(self):
    return self._bounds.size[self.resolution.name_list] / spatial(self.resolution)

prop extrapolation

Expand source code

@property
def extrapolation(self):
    return self._extrapolation

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    raise NotImplementedError

prop face_centers : phiml.math._tensors.Tensor

Expand source code

@property
def face_centers(self) -> Tensor:
    return self._faces.center

prop face_normals : phiml.math._tensors.Tensor

Normal vectors of cell faces, including boundary faces. Shape (elements, ~, vector). For meshes, The vectors point out of the primal cells and into the dual cells.

Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.

Expand source code

@property
def face_normals(self) -> Tensor:
    return self._faces.normal

prop face_shape : phiml.math._shape.Shape

Returns

Expand source code

@property
def face_shape(self) -> Shape:
    return non_channel(self._faces.center)

prop faces : Geometry

Expand source code

@property
def faces(self) -> 'Geometry':
    raise NotImplementedError

prop fill_below

Expand source code

@property
def fill_below(self):
    return self._fill_below

prop grid_bounds

Expand source code

@property
def grid_bounds(self):
    return self._bounds[self.resolution.name_list]

prop height

Expand source code

@property
def height(self):
    return self._height

prop max_dist

Expand source code

@property
def max_dist(self):
    return self._max_dist

prop resolution

Expand source code

@property
def resolution(self):
    return spatial(self._height) - 1

prop shape : phiml.math._shape.Shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self) -> Shape:
    return (self._height.shape - 1) & channel(self._bounds)

prop up

Expand source code

@property
def up(self):
    dims = self._bounds.vector.item_names
    height_unit = vec(**{d: 1 if d == self._hdim else 0 for d in dims})
    return math.where(self._fill_below, height_unit, -height_unit)

prop vertices

Expand source code

@property
def vertices(self):
    hdim = self._hdim
    space = self.vector.item_names
    pos = self.grid_bounds.local_to_global(math.meshgrid(spatial(self._height)) / self.resolution)
    vert = stack({dim: self.height if dim == hdim else pos[dim] for dim in space}, channel('vector'))
    return vert

prop volume : phiml.math._tensors.Tensor

phi.math.Tensor representing the volume of each element. The result retains batch, spatial and instance dimensions.

Expand source code

@property
def volume(self) -> Tensor:
    return math.prod(self.bounding_half_extent() * 2, channel)

Methods

Find the closest surface face of this geometry given a point that can be outside or inside the geometry.

Args

location: Tensor with a single channel dimension called vector. Can have arbitrary other dimensions.

Returns

signed_distance: Scalar signed distance from location to the closest point on the surface. Positive values indicate the point lies outside the geometry, negative values indicate the point lies inside the geometry.
delta: Vector-valued distance vector from location to the closest point on the surface.
normal: Closest surface normal vector.
offset: Min distance of a surface-tangential plane from 0 as a scalar.
face_index: (Optional) An index vector pointing at the closest face.

def approximate_signed_distance(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Computes the approximate distance from location to the surface of the geometry. Locations outside return positive values, inside negative values and zero exactly at the boundary.

When dealing with unions or collections of geometries (instance dimensions), the shortest distance to any instance is returned. This also holds for negative distances.

Args

location: Tensor with one channel dim vector matching the geometry's vector dim.

Returns

Float Tensor

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_half_extent(self) ‑> phiml.math._tensors.Tensor

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def bounding_radius(self) ‑> phiml.math._tensors.Tensor

Returns the radius of a Sphere object that fully encloses this geometry. The sphere is centered at the center of this geometry.

def lies_inside(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def rotated(self, angle: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def sample_uniform(self, *shape: phiml.math._shape.Shape) ‑> phiml.math._tensors.Tensor

Samples uniformly distributed random points inside this volume.

Args

*shape: How many points to sample per individual geometry.

Returns

Tensor containing all dimensions from Geometry.shape, shape as well as a channel dimension vector matching the dimensionality of this Geometry.

def scaled(self, factor: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Scales each individual geometry by factor. The individual center points act as pivots for the operation.

Args

factor: Returns:

def shallow_equals(self, other)

Quick equality check. May return False even if other == self. However, if True is returned, the geometries are guaranteed to be equal.

The shallow_equals() check does not compare all tensor elements but merely checks whether the same tensors are referenced.

class Mesh (vertices: phi.geom._geom.Geometry, elements: phiml.math._tensors.Tensor, element_rank: int, boundaries: Dict[str, Dict[str, slice]], periodic: Sequence[str], face_format: str = 'csc', max_cell_walk: int = None, variable_attrs: Tuple[str, ...] = ('vertices',), value_attrs: Tuple[str, ...] = ())

Unstructured mesh, consisting of vertices and elements.

Use mesh() or mesh_from_numpy() to construct a mesh manually or load_su2() to load one from a file.

Expand source code

@dataclass(frozen=True)
class Mesh(Geometry):
    """
    Unstructured mesh, consisting of vertices and elements.
    
    Use `phi.geom.mesh()` or `phi.geom.mesh_from_numpy()` to construct a mesh manually or `phi.geom.load_su2()` to load one from a file.
    """

    vertices: Geometry
    """ Vertices are represented by a `Geometry` instance with an instance dim. """
    elements: Tensor
    """ elements: Sparse `Tensor` listing ordered vertex indices per element (solid or surface element, depending on `element_rank`).
    Must have one instance dim listing the elements and the corresponding dual dim to `vertices`.
    The vertex count of an element is equal to the number of elements in that row (i.e. summing the dual dim). """
    element_rank: int
    """The spatial rank of the elements. Solid elements have the same as the ambient space, faces one less."""
    boundaries: Dict[str, Dict[str, slice]]
    """Slices to retrieve boundary face values."""
    periodic: Sequence[str]
    """List of axis names that are periodic. Periodic boundaries must be named as axis- and axis+. For example `['x']` will connect the boundaries x- and x+."""
    face_format: str = 'csc'
    """Sparse matrix format for storing quantities that depend on a pair of neighboring elements, e.g. `face_area`, `face_normal`, `face_center`."""
    max_cell_walk: int = None
    """ Maximum number of steps to walk along the element connectivity in order to find a cell, e.g. for sampling at an arbitrary point."""

    variable_attrs: Tuple[str, ...] = ('vertices',)  # PhiML keyword
    value_attrs: Tuple[str, ...] = ()  # PhiML keyword

    def __post_init__(self):
        if spatial(self.elements):
            assert self.elements.dtype.kind == int, f"elements listing vertices must be integer lists but got dtype {self.elements.dtype}"
        else:
            assert self.elements.dtype.kind == bool, f"element matrices must be of type bool but got {self.elements.dtype}"

    @cached_property
    def shape(self) -> Shape:
        return non_dual(self.elements) & channel(self.vertices) & batch(self.vertices)

    @cached_property
    def cell_count(self):
        return instance(self.elements).size

    @cached_property
    def center(self) -> Tensor:
        if self.element_rank == self.spatial_rank:  # Compute volumetric center from faces
            return sum_(self.face_centers * self.face_areas, dual) / sum_(self.face_areas, dual)
        else:  # approximate center from vertices
            return self._vertex_mean

    @cached_property
    def _vertex_mean(self):
        """Mean vertex location per element."""
        vertex_count = sum_(self.elements, instance(self.vertices).as_dual())
        return (self.elements @ self.vertices.center) / vertex_count

    @cached_property
    def face_centers(self) -> Tensor:
        return self._faces['center']

    @property
    def face_areas(self) -> Tensor:
        return self._faces['area']

    @cached_property
    def face_normals(self) -> Tensor:
        if self.element_rank == self.spatial_rank:  # this cannot depend on element centers because that depends on the normals.
            normals = self._faces['normal']
            face_centers = self._faces['center']
            normals_out = normals.vector * (face_centers - self._vertex_mean).vector > 0
            normals = where(normals_out, normals, -normals)
            return normals
        raise NotImplementedError

    @cached_property
    def _faces(self) -> Dict[str, Any]:
        centers, normals, areas, boundary_slices = build_faces(self.vertices.center, self.elements, self.boundaries, self.element_rank, self.periodic, self._vertex_mean, self.face_format)
        return {
            'center': centers,
            'normal': normals,
            'area': areas,
            'boundary_slices': boundary_slices,
        }

    @property
    def face_shape(self) -> Shape:
        return instance(self.elements) & dual

    @property
    def sets(self):
        return {
            'center': non_batch(self)-'vector',
            'vertex': instance(self.vertices),
            '~vertex': dual(self.elements)
        }

    def get_points(self, set_key: str) -> Tensor:
        if set_key == 'vertex':
            return self.vertices.center
        elif set_key == '~vertex':
            return si2d(self.vertices.center)
        else:
            return Geometry.get_points(self, set_key)

    def get_boundary(self, set_key: str) -> Dict[str, Dict[str, slice]]:
        if set_key in ['vertex', '~vertex']:
            return {}
        return Geometry.get_boundary(self, set_key)

    @property
    def boundary_elements(self) -> Dict[str, Dict[str, slice]]:
        return {}

    @property
    def boundary_faces(self) -> Dict[str, Dict[str, slice]]:
        return self._faces['boundary_slices']

    @property
    def all_boundary_faces(self) -> Dict[str, slice]:
        return {self.face_shape.dual.name: slice(instance(self).volume, None)}
    
    @property
    def interior_faces(self) -> Dict[str, slice]:
        return {self.face_shape.dual.name: slice(0, instance(self).volume)}

    def pad_boundary(self, value: Tensor, widths: Dict[str, Dict[str, slice]] = None, mode: Extrapolation or Tensor or Number = 0, **kwargs) -> Tensor:
        mode = as_extrapolation(mode)
        if self.face_shape.dual.name not in value.shape:
            value = rename_dims(value, instance, self.face_shape.dual)
        else:
            raise NotImplementedError
        if widths is None:
            widths = self.boundary_faces
        if isinstance(widths, (tuple, list)):
            if len(widths) == 0 or isinstance(widths[0], dict):  # add sliced-off slices
                pass
        dim = next(iter(next(iter(widths.values()))))
        slices = [slice(0, value.shape.get_size(dim))]
        values = [value]
        connectivity = self.connectivity
        for name, b_slice in widths.items():
            if b_slice[dim].stop - b_slice[dim].start > 0:
                slices.append(b_slice[dim])
                values.append(mode.sparse_pad_values(value, connectivity[b_slice], name, mesh=self, **kwargs))
        perm = np.argsort([s.start for s in slices])
        ordered_pieces = [values[i] for i in perm]
        return concat(ordered_pieces, dim, expand_values=True)

    @cached_property
    def cell_connectivity(self) -> Tensor:
        """
        Returns a bool-like matrix whose non-zero entries denote connected elements.
        In meshes or grids, elements are connected if they share a face in 3D, an edge in 2D, or a vertex in 1D.

        Returns:
            `Tensor` of shape (elements, ~elements)
        """
        return self.connectivity[self.interior_faces]

    @cached_property
    def boundary_connectivity(self) -> Tensor:
        return self.connectivity[self.all_boundary_faces]

    @cached_property
    def distance_matrix(self):
        return vec_length(pairwise_distances(self.center, edges=self.cell_connectivity, format='as edges', default=None))

    def faces_to_vertices(self, values: Tensor, reduce=sum):
        v = stored_values(values, invalid='keep')  # ToDo replace this once PhiML has support for dense instance dims and sparse scatter
        i = stored_values(self.face_vertices, invalid='keep')
        i = rename_dims(i, channel, instance)
        out_shape = non_channel(self.vertices) & shape(values).without(self.face_shape)
        return scatter(out_shape, i, v, mode=reduce, outside_handling='undefined')

    @cached_property
    def _cell_deltas(self):
        bounds = bounding_box(self.vertices)
        periodic = {dim[:-len('[::-1]')] if dim.endswith('[::-1]') else dim: dim.endswith('[::-1]') for dim in self.periodic}
        is_periodic = dim_mask(self.vector.item_names, tuple(periodic))
        return pairwise_distances(self.center, format=self.cell_connectivity, periodic=is_periodic, domain=(bounds.lower, bounds.upper))

    @cached_property
    def relative_face_distance(self):
        """|face_center - center| / |neighbor_center - center|"""
        cell_distances = vec_length(self._cell_deltas)
        assert (cell_distances > 0).all, f"All cells must have distance > 0 but found 0 distance at {nonzero(cell_distances == 0)}"
        face_distances = vec_length(self.face_centers[self.interior_faces] - self.center)
        return concat([face_distances / cell_distances, self.boundary_connectivity], self.face_shape.dual)

    @cached_property
    def neighbor_offsets(self):
        """Returns shift vector to neighbor centroids and boundary faces."""
        boundary_deltas = (self.face_centers - self.center)[self.all_boundary_faces]
        assert (vec_length(boundary_deltas) > 0).all, f"All boundary faces must be separated from the cell centers but 0 distance at the following {channel(stored_indices(boundary_deltas)).item_names[0]}:\n{nonzero(vec_length(boundary_deltas) == 0):full}"
        return concat([self._cell_deltas, boundary_deltas], self.face_shape.dual)

    @cached_property
    def neighbor_distances(self):
        return vec_length(self.neighbor_offsets)

    @property
    def faces(self) -> 'Geometry':
        """
        Assembles information about the boundaries of the elements that make up the surface.
        For 2D elements, the faces are edges, for 3D elements, the faces are planar elements.

        Returns:
            center: Center of face connecting a pair of elements. Shape (~elements, elements, vector).
                Returns 0-vectors for unconnected elements.
            area: Area of face connecting a pair of elements. Shape (~elements, elements).
                Returns 0 for unconnected elements.
            normal: Normal vector of face connecting a pair of elements. Shape (~elements, elements, vector).
                Unconnected elements are assigned the vector 0.
                The vector points out of polygon and into ~polygon.
        """
        return Point(self.face_centers)

    @property
    def connectivity(self) -> Tensor:
        return self.element_connectivity

    @cached_property
    def element_connectivity(self) -> Tensor:
        if self.element_rank == self.spatial_rank:
            if is_sparse(self.face_areas):
                return tensor_like(self.face_areas, True)
            else:
                return self.face_areas > 0
        else:  # fallback with no boundaries
            coo = to_format(self.elements, 'coo').numpy()
            connected_elements = coo @ coo.T
            connected_elements.data = np.ones_like(connected_elements.data)
            element_connectivity = wrap(connected_elements, instance(self.elements), instance(self.elements).as_dual())
            return element_connectivity

    @cached_property
    def vertex_connectivity(self) -> Tensor:
        if isinstance(self.vertices, Graph):
            return self.vertices.connectivity
        elif self.element_rank <= 2:
            def single_vertex_connectivity(elements: Tensor):
                indices = stored_indices(elements).index[dual(elements).name]
                idx1 = indices.numpy()
                v_count = sum_(elements, dual).numpy()
                ptr_end = np.cumsum(v_count)
                roll = np.arange(idx1.size) + 1
                roll[ptr_end-1] = ptr_end - v_count
                idx2 = idx1[roll]
                v_conn = coo_matrix((np.ones(idx1.size, dtype=bool), (idx1, idx2)), shape=(dual(elements).size,)*2).tocsr()
                return wrap(v_conn, dual(elements).as_instance(), dual(elements))
            return math.map(single_vertex_connectivity, self.elements, dims=batch)
        raise NotImplementedError

    @cached_property
    def vertex_graph(self) -> Graph:
        return self.vertices if isinstance(self.vertices, Graph) else graph(self.vertices, self.vertex_connectivity)

    def filter_unused_vertices(self) -> 'Mesh':
        coo = to_format(self.elements, 'coo').numpy()
        has_element = np.asarray(coo.sum(0) > 0)[0]
        new_index = np.cumsum(has_element) - 1
        new_index_t = wrap(new_index, dual(self.elements))
        has_element = wrap(has_element, instance(self.vertices))
        has_element_d = si2d(has_element)
        vertices = self.vertices[has_element]
        v_normals = self.vertex_normals[has_element_d]
        vertex_connectivity = None
        # if self._vertex_connectivity is not None:
        #     vertex_connectivity = stored_indices(self._vertex_connectivity).index.as_batch()
        #     vertex_connectivity = new_index_t[{dual: vertex_connectivity}].index.as_channel()
        #     vertex_connectivity = sparse_tensor(vertex_connectivity, stored_values(self._vertex_connectivity), non_batch(self._vertex_connectivity).with_sizes(instance(vertices).size), False)
        if isinstance(self.elements, CompactSparseTensor):
            indices = new_index_t[{dual: self.elements._indices}]
            elements = CompactSparseTensor(indices, self.elements._values, self.elements._compressed_dims.with_size(instance(vertices).volume), self.elements._indices_constant, self.elements._matrix_rank)
        else:
            filtered_coo = coo_matrix((coo.data, (coo.row, new_index)), shape=(instance(self.elements).volume, instance(vertices).volume))  # ToDo keep sparse format
            elements = wrap(filtered_coo, self.elements.shape.without_sizes())
        return Mesh(vertices, elements, self.element_rank, self.boundaries, self.periodic, self.face_format, self.max_cell_walk, self.variable_attrs, self.value_attrs)

    @property
    def volume(self) -> Tensor:
        if self.element_rank == 2:
            if instance(self.elements).volume > 0:
                three_vertices = nonzero(self.elements, 3, list_dims=dual)
                v1, v2, v3 = unstack(self.vertices.center[{instance: three_vertices}], dual)
                cross_area = vec_length(cross(v2-v1, v3-v1))
                vertex_count = math.sum(self.elements, dual)
                fac = where(vertex_count == 3, 0.5, 1)  # tri, quad, ...
                return fac * cross_area
            else:
                return zeros(instance(self.vertices))  # empty mesh
        elif self.element_rank == self.spatial_rank:
            vol_contributions = (self.face_centers.vector @ self.face_normals.vector) * self.face_areas
            return sum_(vol_contributions, dual) / self.spatial_rank
        raise NotImplementedError


    @property
    def normals(self) -> Tensor:
        """Extrinsic element normal space. This is a 0D vector for solid elements and 1D for surface elements."""
        if self.element_rank == 2:
            three_vertices = nonzero(self.elements, 3, list_dims=dual)
            v1, v2, v3 = unstack(self.vertices.center[{instance: three_vertices}], dual)
            return vec_normalize(cross(v2 - v1, v3 - v1))
        raise NotImplementedError

    @property
    def vertex_normals(self) -> Tensor:
        v_normals = mean(self.elements * self.normals, instance)  # (~vertices,vector)
        return vec_normalize(v_normals)

    @property
    def vertex_positions(self) -> Tensor:
        """Lists the vertex centers along the corresponding dual dim to `self.vertices.center`."""
        return si2d(self.vertices.center)

    def lies_inside(self, location: Tensor) -> Tensor:
        idx = find_closest(self.center, location)
        for i in range(self.max_cell_walk):
            idx, leaves_mesh, is_outside, *_ = self.cell_walk_towards(location, idx, allow_exit=i == self.max_cell_walk - 1)
        return ~(leaves_mesh & is_outside)

    def approximate_signed_distance(self, location: Union[Tensor, tuple]) -> Tensor:
        if self.element_rank == 2 and self.spatial_rank == 3:
            closest_elem = find_closest(self.center, location)
            center = self.center[closest_elem]
            normal = self.normals[closest_elem]
            return plane_sgn_dist(center, normal, location)
        idx = find_closest(self.center, location)
        for i in range(self.max_cell_walk):
            idx, leaves_mesh, is_outside, distances, nb_idx = self.cell_walk_towards(location, idx, allow_exit=False)
        return math.max(distances, dual)

    def approximate_closest_surface(self, location: Tensor) -> Tuple[Tensor, Tensor, Tensor, Tensor, Tensor]:
        if self.element_rank == 2 and self.spatial_rank == 3:
            closest_elem = find_closest(self.center, location)
            center = self.center[closest_elem]
            normal = self.normals[closest_elem]
            face_size = sqrt(self.volume) * 4
            size = face_size[closest_elem]
            sgn_dist = plane_sgn_dist(center, normal, location)
            delta = center - location  # this is not accurate...
            outward = where(abs(sgn_dist) < size, normal, vec_normalize(delta))
            return sgn_dist, delta, outward, None, closest_elem
        # idx = find_closest(self.center, location)
        # for i in range(self.max_cell_walk):
        #     idx, leaves_mesh, is_outside, distances, nb_idx = self.cell_walk_towards(location, idx, allow_exit=False)
        # sgn_dist = max(distances, dual)
        # cell_normals = self.face_normals[idx]
        # normal = cell_normals[{dual: nb_idx}]
        # return sgn_dist, delta, normal, offset, face_index
        raise NotImplementedError

    def cell_walk_towards(self, location: Tensor, start_cell_idx: Tensor, allow_exit=False):
        """
        If `location` is not within the cell at index `from_cell_idx`, moves to a closer neighbor cell.

        Args:
            location: Target location as `Tensor`.
            start_cell_idx: Index of starting cell. Must be a valid cell index.
            allow_exit: If `True`, returns an invalid index for points outside the mesh, otherwise keeps the current index.

        Returns:
            index: Index of the neighbor cell or starting cell.
            leaves_mesh: Whether the walk crossed the mesh boundary. Then `index` is invalid. This is only possible if `allow_exit` is true.
            is_outside: Whether `location` was outside the cell at index `start_cell_idx`.
        """
        closest_normals = self.face_normals[start_cell_idx]
        closest_face_centers = self.face_centers[start_cell_idx]
        offsets = closest_normals.vector @ closest_face_centers.vector  # this dot product could be cashed in the mesh
        distances = closest_normals.vector @ location.vector - offsets
        is_outside = math.any(distances > 0, dual)
        nb_idx = argmax(distances, dual).index[0]  # cell index or boundary face index
        leaves_mesh = nb_idx >= instance(self).volume
        next_idx = where(is_outside & (~leaves_mesh | allow_exit), nb_idx, start_cell_idx)
        return next_idx, leaves_mesh, is_outside, distances, nb_idx

    def sample_uniform(self, *shape: Shape) -> Tensor:
        raise NotImplementedError

    def bounding_radius(self) -> Tensor:
        center = self.elements * self.center
        vert_pos = rename_dims(self.vertices.center, instance, dual)
        dist_to_vert = vec_length(vert_pos - center)
        max_dist = math.max(dist_to_vert, dual)
        return max_dist

    def bounding_half_extent(self) -> Tensor:
        center = self.elements * self.center
        vert_pos = rename_dims(self.vertices.center, instance, dual)
        max_delta = math.max(abs(vert_pos - center), dual)
        return max_delta

    def bounding_box(self) -> 'BaseBox':
        return self.vertices.bounding_box()

    @property
    def bounds(self):
        return Box(math.min(self.vertices.center, instance), math.max(self.vertices.center, instance))

    def at(self, center: Tensor) -> 'Mesh':
        if instance(self.elements) in center.shape:
            raise NotImplementedError("Setting Mesh positions only supported for vertices, not elements")
        if dual(self.elements) in center.shape:
            center = rename_dims(center, dual, instance(self.vertices))
        if instance(self.vertices) in center.shape:
            vertices = self.vertices.at(center)
            return Mesh(vertices, self.elements, self.element_rank, self.boundaries, self.periodic, self.face_format, self.max_cell_walk, self.variable_attrs, self.value_attrs)
        else:
            return self.shifted(center - self.bounds.center)

    def shifted(self, delta: Tensor) -> 'Mesh':
        if instance(self.elements) in delta.shape:
            raise NotImplementedError("Shifting Mesh positions only supported for vertices, not elements")
        if dual(self.elements) in delta.shape:
            delta = rename_dims(delta, dual, instance(self.vertices))
        if instance(self.vertices) in delta.shape:
            vertices = self.vertices.shifted(delta)
            return Mesh(vertices, self.elements, self.element_rank, self.boundaries, self.periodic, self.face_format, self.max_cell_walk, self.variable_attrs, self.value_attrs)
        else:  # shift everything
            # ToDo transfer cached properties
            # copy: center+delta, normals, volume, face_centers+delta, face_areas, face_normals, vertex_normals, vertex_connectivity, element_connectivity
            vertices = self.vertices.shifted(delta)
            return Mesh(vertices, self.elements, self.element_rank, self.boundaries, self.periodic, self.face_format, self.max_cell_walk, self.variable_attrs, self.value_attrs)

    def rotated(self, angle: Union[float, Tensor]) -> 'Geometry':
        raise NotImplementedError

    def scaled(self, factor: float | Tensor) -> 'Geometry':
        pivot = self.bounds.center
        vertices = scale(self.vertices, factor, pivot)
        center = scale(Point(self.center), factor, pivot).center
        volume = self.volume * factor**self.element_rank if self.volume is not None else None
        face_areas = None
        return Mesh(vertices, self.elements, self.element_rank, self.boundaries, center, volume, self.normals, self.face_centers, self.face_normals, face_areas, self.face_vertices, self._vertex_normals, self._vertex_connectivity, self._element_connectivity, self.max_cell_walk)

    def __getitem__(self, item):
        item: dict = slicing_dict(self, item)
        assert not spatial(self.elements).only(tuple(item)), f"Cannot slice vertex lists ('{spatial(self.elements)}') but got slicing dict {item}"
        assert not instance(self.vertices).only(tuple(item)), f"Slicing by vertex indices ('{instance(self.vertices)}') not supported but got slicing dict {item}"
        return getitem_dataclass(self, item, keepdims=[self.shape.instance.name, 'vector'])

    def __repr__(self):
        return Geometry.__repr__(self)

Ancestors

phi.geom._geom.Geometry

Class variables

var boundaries : Dict[str, Dict[str, slice]]: Slices to retrieve boundary face values.
var element_rank : int: The spatial rank of the elements. Solid elements have the same as the ambient space, faces one less.
var elements : phiml.math._tensors.Tensor: elements: Sparse Tensor listing ordered vertex indices per element (solid or surface element, depending on element_rank). Must have one instance dim listing the elements and the corresponding dual dim to vertices. The vertex count of an element is equal to the number of elements in that row (i.e. summing the dual dim).
var face_format : str: Sparse matrix format for storing quantities that depend on a pair of neighboring elements, e.g. face_area, face_normal, face_center.
var max_cell_walk : int: Maximum number of steps to walk along the element connectivity in order to find a cell, e.g. for sampling at an arbitrary point.
var periodic : Sequence[str]: List of axis names that are periodic. Periodic boundaries must be named as axis- and axis+. For example ['x'] will connect the boundaries x- and x+.
var value_attrs : Tuple[str, ...]
var variable_attrs : Tuple[str, ...]
var vertices : phi.geom._geom.Geometry: Vertices are represented by a Geometry instance with an instance dim.

Instance variables

prop all_boundary_faces : Dict[str, slice]

Expand source code

@property
def all_boundary_faces(self) -> Dict[str, slice]:
    return {self.face_shape.dual.name: slice(instance(self).volume, None)}

var boundary_connectivity

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

prop boundary_elements : Dict[str, Dict[str, slice]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[str, Dict[str, slice]]:
    return {}

prop boundary_faces : Dict[str, Dict[str, slice]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[str, Dict[str, slice]]:
    return self._faces['boundary_slices']

prop bounds

Expand source code

@property
def bounds(self):
    return Box(math.min(self.vertices.center, instance), math.max(self.vertices.center, instance))

var cell_connectivity

Returns a bool-like matrix whose non-zero entries denote connected elements. In meshes or grids, elements are connected if they share a face in 3D, an edge in 2D, or a vertex in 1D.

Returns

Tensor of shape (elements, ~elements)

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

var cell_count

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

var center

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

prop connectivity : phiml.math._tensors.Tensor

Expand source code

@property
def connectivity(self) -> Tensor:
    return self.element_connectivity

var distance_matrix

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

var element_connectivity

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    return self._faces['area']

var face_centers

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

var face_normals

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

prop face_shape : phiml.math._shape.Shape

Returns

Expand source code

@property
def face_shape(self) -> Shape:
    return instance(self.elements) & dual

prop faces : Geometry

Assembles information about the boundaries of the elements that make up the surface. For 2D elements, the faces are edges, for 3D elements, the faces are planar elements.

Returns

center: Center of face connecting a pair of elements. Shape (~elements, elements, vector). Returns 0-vectors for unconnected elements.
area: Area of face connecting a pair of elements. Shape (~elements, elements). Returns 0 for unconnected elements.
normal: Normal vector of face connecting a pair of elements. Shape (~elements, elements, vector). Unconnected elements are assigned the vector 0. The vector points out of polygon and into ~polygon.

Expand source code

@property
def faces(self) -> 'Geometry':
    """
    Assembles information about the boundaries of the elements that make up the surface.
    For 2D elements, the faces are edges, for 3D elements, the faces are planar elements.

    Returns:
        center: Center of face connecting a pair of elements. Shape (~elements, elements, vector).
            Returns 0-vectors for unconnected elements.
        area: Area of face connecting a pair of elements. Shape (~elements, elements).
            Returns 0 for unconnected elements.
        normal: Normal vector of face connecting a pair of elements. Shape (~elements, elements, vector).
            Unconnected elements are assigned the vector 0.
            The vector points out of polygon and into ~polygon.
    """
    return Point(self.face_centers)

prop interior_faces : Dict[str, slice]

Expand source code

@property
def interior_faces(self) -> Dict[str, slice]:
    return {self.face_shape.dual.name: slice(0, instance(self).volume)}

var neighbor_distances

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

var neighbor_offsets

Returns shift vector to neighbor centroids and boundary faces.

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

prop normals : phiml.math._tensors.Tensor

Extrinsic element normal space. This is a 0D vector for solid elements and 1D for surface elements.

Expand source code

@property
def normals(self) -> Tensor:
    """Extrinsic element normal space. This is a 0D vector for solid elements and 1D for surface elements."""
    if self.element_rank == 2:
        three_vertices = nonzero(self.elements, 3, list_dims=dual)
        v1, v2, v3 = unstack(self.vertices.center[{instance: three_vertices}], dual)
        return vec_normalize(cross(v2 - v1, v3 - v1))
    raise NotImplementedError

var relative_face_distance

|face_center - center| / |neighbor_center - center|

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

prop sets

Expand source code

@property
def sets(self):
    return {
        'center': non_batch(self)-'vector',
        'vertex': instance(self.vertices),
        '~vertex': dual(self.elements)
    }

var shape

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

var vertex_connectivity

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

var vertex_graph

Expand source code

def __get__(self, instance, owner=None):
    if instance is None:
        return self
    if self.attrname is None:
        raise TypeError(
            "Cannot use cached_property instance without calling __set_name__ on it.")
    try:
        cache = instance.__dict__
    except AttributeError:  # not all objects have __dict__ (e.g. class defines slots)
        msg = (
            f"No '__dict__' attribute on {type(instance).__name__!r} "
            f"instance to cache {self.attrname!r} property."
        )
        raise TypeError(msg) from None
    val = cache.get(self.attrname, _NOT_FOUND)
    if val is _NOT_FOUND:
        val = self.func(instance)
        try:
            cache[self.attrname] = val
        except TypeError:
            msg = (
                f"The '__dict__' attribute on {type(instance).__name__!r} instance "
                f"does not support item assignment for caching {self.attrname!r} property."
            )
            raise TypeError(msg) from None
    return val

prop vertex_normals : phiml.math._tensors.Tensor

Expand source code

@property
def vertex_normals(self) -> Tensor:
    v_normals = mean(self.elements * self.normals, instance)  # (~vertices,vector)
    return vec_normalize(v_normals)

prop vertex_positions : phiml.math._tensors.Tensor

Lists the vertex centers along the corresponding dual dim to self.vertices.center.

Expand source code

@property
def vertex_positions(self) -> Tensor:
    """Lists the vertex centers along the corresponding dual dim to `self.vertices.center`."""
    return si2d(self.vertices.center)

prop volume : phiml.math._tensors.Tensor

phi.math.Tensor representing the volume of each element. The result retains batch, spatial and instance dimensions.

Expand source code

@property
def volume(self) -> Tensor:
    if self.element_rank == 2:
        if instance(self.elements).volume > 0:
            three_vertices = nonzero(self.elements, 3, list_dims=dual)
            v1, v2, v3 = unstack(self.vertices.center[{instance: three_vertices}], dual)
            cross_area = vec_length(cross(v2-v1, v3-v1))
            vertex_count = math.sum(self.elements, dual)
            fac = where(vertex_count == 3, 0.5, 1)  # tri, quad, ...
            return fac * cross_area
        else:
            return zeros(instance(self.vertices))  # empty mesh
    elif self.element_rank == self.spatial_rank:
        vol_contributions = (self.face_centers.vector @ self.face_normals.vector) * self.face_areas
        return sum_(vol_contributions, dual) / self.spatial_rank
    raise NotImplementedError

Methods

Find the closest surface face of this geometry given a point that can be outside or inside the geometry.

Args

location: Tensor with a single channel dimension called vector. Can have arbitrary other dimensions.

Returns

signed_distance: Scalar signed distance from location to the closest point on the surface. Positive values indicate the point lies outside the geometry, negative values indicate the point lies inside the geometry.
delta: Vector-valued distance vector from location to the closest point on the surface.
normal: Closest surface normal vector.
offset: Min distance of a surface-tangential plane from 0 as a scalar.
face_index: (Optional) An index vector pointing at the closest face.

def approximate_signed_distance(self, location: phiml.math._tensors.Tensor | tuple) ‑> phiml.math._tensors.Tensor

Computes the approximate distance from location to the surface of the geometry. Locations outside return positive values, inside negative values and zero exactly at the boundary.

When dealing with unions or collections of geometries (instance dimensions), the shortest distance to any instance is returned. This also holds for negative distances.

Args

location: Tensor with one channel dim vector matching the geometry's vector dim.

Returns

Float Tensor

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._mesh.Mesh

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_box(self) ‑> phi.geom._box.BaseBox

Returns the approximately smallest axis-aligned box that contains this Geometry. The center of the box may not be equal to self.center.

Returns

Box or Cuboid that fully contains this Geometry.

def bounding_half_extent(self) ‑> phiml.math._tensors.Tensor

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def bounding_radius(self) ‑> phiml.math._tensors.Tensor

Returns the radius of a Sphere object that fully encloses this geometry. The sphere is centered at the center of this geometry.

def cell_walk_towards(self, location: phiml.math._tensors.Tensor, start_cell_idx: phiml.math._tensors.Tensor, allow_exit=False)

If location is not within the cell at index from_cell_idx, moves to a closer neighbor cell.

Args

location: Target location as Tensor.
start_cell_idx: Index of starting cell. Must be a valid cell index.
allow_exit: If True, returns an invalid index for points outside the mesh, otherwise keeps the current index.

Returns

index: Index of the neighbor cell or starting cell.
leaves_mesh: Whether the walk crossed the mesh boundary. Then index is invalid. This is only possible if allow_exit is true.
is_outside: Whether location was outside the cell at index start_cell_idx.

def faces_to_vertices(self, values: phiml.math._tensors.Tensor, reduce=<built-in function sum>)

def filter_unused_vertices(self) ‑> phi.geom._mesh.Mesh

def get_boundary(self, set_key: str) ‑> Dict[str, Dict[str, slice]]

def get_points(self, set_key: str) ‑> phiml.math._tensors.Tensor

def lies_inside(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def pad_boundary(self, value: phiml.math._tensors.Tensor, widths: Dict[str, Dict[str, slice]] = None, mode: phiml.math.extrapolation.Extrapolation = 0, **kwargs) ‑> phiml.math._tensors.Tensor

def rotated(self, angle: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def sample_uniform(self, *shape: phiml.math._shape.Shape) ‑> phiml.math._tensors.Tensor

Samples uniformly distributed random points inside this volume.

Args

*shape: How many points to sample per individual geometry.

Returns

Tensor containing all dimensions from Geometry.shape, shape as well as a channel dimension vector matching the dimensionality of this Geometry.

def scaled(self, factor: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Scales each individual geometry by factor. The individual center points act as pivots for the operation.

Args

factor: Returns:

def shifted(self, delta: phiml.math._tensors.Tensor) ‑> phi.geom._mesh.Mesh

Returns a translated version of this geometry.

See Also: Geometry.at().

Args

delta: direction vector
delta: Tensor:

Returns

Geometry: shifted geometry

class Point (location: phiml.math._tensors.Tensor)

Points have zero volume and are determined by a single location. An instance of Point represents a single n-dimensional point or a batch of points.

Expand source code

class Point(Geometry):
    """
    Points have zero volume and are determined by a single location.
    An instance of `Point` represents a single n-dimensional point or a batch of points.
    """

    def __init__(self, location: math.Tensor):
        assert 'vector' in location.shape, "location must have a vector dimension"
        assert location.shape.get_item_names('vector') is not None, "Vector dimension needs to list spatial dimension as item names."
        self._location = location
        self._shape = self._location.shape

    def __variable_attrs__(self):
        return '_location',

    def __value_attrs__(self):
        return '_location',

    def __with_attrs__(self, **updates):
        if '_location' in updates:
            result = Point.__new__(Point)
            result._location = updates['_location']
            result._shape = result._location.shape if result._location is not None else self._shape
            return result
        else:
            return self

    @property
    def center(self) -> Tensor:
        return self._location

    @property
    def shape(self) -> Shape:
        return self._shape

    @property
    def faces(self) -> 'Geometry':
        return self

    def unstack(self, dimension: str) -> tuple:
        return tuple(Point(loc) for loc in math.unstack(self._location, dimension))

    def lies_inside(self, location: Tensor) -> Tensor:
        return expand(math.wrap(False), shape(location).without('vector'))

    def approximate_signed_distance(self, location: Union[Tensor, tuple]) -> Tensor:
        return math.vec_abs(location - self._location)

    def bounding_radius(self) -> Tensor:
        return math.zeros()

    def bounding_half_extent(self) -> Tensor:
        return tensor_like(self.center, 0)

    def at(self, center: Tensor) -> 'Geometry':
        return Point(center)

    def rotated(self, angle) -> 'Geometry':
        return self

    @property
    def volume(self) -> Tensor:
        return math.wrap(0)

    def sample_uniform(self, *shape: math.Shape) -> Tensor:
        raise NotImplementedError

    def scaled(self, factor: Union[float, Tensor]) -> 'Geometry':
        return self

    @property
    def face_centers(self) -> Tensor:
        return self._location

    @property
    def face_areas(self) -> Tensor:
        return expand(0, self.face_shape)

    @property
    def face_normals(self) -> Tensor:
        raise AssertionError(f"Points have no normals")

    @property
    def boundary_elements(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
        return {}

    @property
    def boundary_faces(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
        return {}

    @property
    def face_shape(self) -> Shape:
        return self.shape

    @property
    def corners(self):
        return self._location

    def __getitem__(self, item):
        return Point(self._location[_keep_vector(slicing_dict(self, item))])

Ancestors

phi.geom._geom.Geometry

Instance variables

prop boundary_elements : Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
    return {}

prop boundary_faces : Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
    return {}

prop center : phiml.math._tensors.Tensor

Center location in single channel dimension.

Expand source code

@property
def center(self) -> Tensor:
    return self._location

prop corners

Returns

Expand source code

@property
def corners(self):
    return self._location

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    return expand(0, self.face_shape)

prop face_centers : phiml.math._tensors.Tensor

Expand source code

@property
def face_centers(self) -> Tensor:
    return self._location

prop face_normals : phiml.math._tensors.Tensor

Normal vectors of cell faces, including boundary faces. Shape (elements, ~, vector). For meshes, The vectors point out of the primal cells and into the dual cells.

Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.

Expand source code

@property
def face_normals(self) -> Tensor:
    raise AssertionError(f"Points have no normals")

prop face_shape : phiml.math._shape.Shape

Returns

Expand source code

@property
def face_shape(self) -> Shape:
    return self.shape

prop faces : Geometry

Expand source code

@property
def faces(self) -> 'Geometry':
    return self

prop shape : phiml.math._shape.Shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self) -> Shape:
    return self._shape

prop volume : phiml.math._tensors.Tensor

phi.math.Tensor representing the volume of each element. The result retains batch, spatial and instance dimensions.

Expand source code

@property
def volume(self) -> Tensor:
    return math.wrap(0)

Methods

def approximate_signed_distance(self, location: phiml.math._tensors.Tensor | tuple) ‑> phiml.math._tensors.Tensor

Computes the approximate distance from location to the surface of the geometry. Locations outside return positive values, inside negative values and zero exactly at the boundary.

When dealing with unions or collections of geometries (instance dimensions), the shortest distance to any instance is returned. This also holds for negative distances.

Args

location: Tensor with one channel dim vector matching the geometry's vector dim.

Returns

Float Tensor

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_half_extent(self) ‑> phiml.math._tensors.Tensor

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def bounding_radius(self) ‑> phiml.math._tensors.Tensor

Returns the radius of a Sphere object that fully encloses this geometry. The sphere is centered at the center of this geometry.

def lies_inside(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def rotated(self, angle) ‑> phi.geom._geom.Geometry

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def sample_uniform(self, *shape: phiml.math._shape.Shape) ‑> phiml.math._tensors.Tensor

Samples uniformly distributed random points inside this volume.

Args

*shape: How many points to sample per individual geometry.

Returns

Tensor containing all dimensions from Geometry.shape, shape as well as a channel dimension vector matching the dimensionality of this Geometry.

def scaled(self, factor: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Scales each individual geometry by factor. The individual center points act as pivots for the operation.

Args

factor: Returns:

def unstack(self, dimension: str) ‑> tuple

Unstacks this Geometry along the given dimension. The shapes of the returned geometries are reduced by dimension.

Args

dimension: dimension along which to unstack

Returns

geometries: tuple of length equal to geometry.shape.get_size(dimension)

class SDF (sdf: Callable, out_shape=None, bounds: phi.geom._box.BaseBox = None, center: phiml.math._tensors.Tensor = None, volume: phiml.math._tensors.Tensor = None, bounding_radius: phiml.math._tensors.Tensor = None, sdf_and_grad: Callable = None)

Function-based signed distance field. Negative values lie inside the geometry, the 0-level represents the surface.

Args

sdf: SDF function. First argument is a phiml.math.Tensor with a vector channel dim.
bounds: Grid limits. The bounds fully enclose all virtual cells.
center: (Optional) Geometry center point. Will be computed otherwise.
volume: (Optional) Geometry volume. Will be computed otherwise.
bounding_radius: (Optional) Geometry bounding radius around center. Will be computed otherwise.

Expand source code

class SDF(Geometry):
    """
    Function-based signed distance field.
    Negative values lie inside the geometry, the 0-level represents the surface.
    """
    def __init__(self, sdf: Callable, out_shape=None, bounds: BaseBox = None, center: Tensor = None, volume: Tensor = None, bounding_radius: Tensor = None, sdf_and_grad: Callable = None):
        """
        Args:
            sdf: SDF function. First argument is a `phiml.math.Tensor` with a `vector` channel dim.
            bounds: Grid limits. The bounds fully enclose all virtual cells.
            center: (Optional) Geometry center point. Will be computed otherwise.
            volume: (Optional) Geometry volume. Will be computed otherwise.
            bounding_radius: (Optional) Geometry bounding radius around center. Will be computed otherwise.
        """
        self._sdf = sdf
        if out_shape is not None:
            self._out_shape = out_shape or math.EMPTY_SHAPE
        else:
            dims = channel([bounds, center, bounding_radius])
            assert 'vector' in dims, f"If out_shape is not specified, either bounds, center or bounding_radius must be given."
            self._out_shape = sdf(math.zeros(dims['vector'])).shape
        self._bounds = bounds
        if sdf_and_grad is not None:
            self._grad = sdf_and_grad
        else:
            self._grad = math.gradient(sdf, wrt=0, get_output=True)
        if center is not None:
            self._center = center
        else:
            self._center = bounds.center
        if volume is not None:
            self._volume = volume
        else:
            self._volume = None
        if bounding_radius is not None:
            self._bounding_radius = bounding_radius
        else:
            self._bounding_radius = self._bounds.bounding_radius()

    def __call__(self, location, *aux_args, **aux_kwargs):
        native_loc = not isinstance(location, Tensor)
        if native_loc:
            location = math.wrap(location, instance('points'), self.shape['vector'])
        sdf_val: Tensor = self._sdf(location, *aux_args, **aux_kwargs)
        return sdf_val.native() if native_loc else sdf_val

    def __variable_attrs__(self):
        return '_bounds', '_center', '_volume', '_bounding_radius'

    def __value_attrs__(self):
        return ()

    @property
    def bounds(self) -> BaseBox:
        return self._bounds

    @property
    def size(self):
        return self._bounds.size

    @property
    def center(self) -> Tensor:
        return self._center

    @property
    def shape(self) -> Shape:
        return self._out_shape & self._bounds.shape

    @property
    def volume(self) -> Tensor:
        return self._volume

    @property
    def faces(self) -> 'Geometry':
        raise NotImplementedError(f"SDF does not support faces")

    @property
    def face_centers(self) -> Tensor:
        raise NotImplementedError(f"SDF does not support faces")

    @property
    def face_areas(self) -> Tensor:
        raise NotImplementedError(f"SDF does not support faces")

    @property
    def face_normals(self) -> Tensor:
        raise NotImplementedError(f"SDF does not support faces")

    @property
    def boundary_elements(self) -> Dict[Any, Dict[str, slice]]:
        return {}

    @property
    def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
        return {}

    @property
    def face_shape(self) -> Shape:
        return math.EMPTY_SHAPE

    @property
    def corners(self) -> Tensor:
        raise NotImplementedError(f"SDF does not support corners")

    def lies_inside(self, location: Tensor) -> Tensor:
        sdf = self._sdf(location)
        return sdf <= 0

    def approximate_closest_surface(self, location: Tensor, refine_iter=0) -> Tuple[Tensor, Tensor, Tensor, Tensor, Tensor]:
        sgn_dist, outward = self._grad(location)
        closest = location - sgn_dist * outward
        if not refine_iter:
            _, normal = self._grad(closest)
        else:
            for i in range(refine_iter):
                sgn_dist, outward = self._grad(closest)
                closest -= sgn_dist * outward
            normal = outward
        offset = None
        face_index = None
        return sgn_dist, closest - location, normal, offset, face_index

    def sdf_and_gradient(self, location: Tensor, refine_iter=0) -> Tuple[Tensor, Tensor]:
        if not refine_iter:
            sgn_dist, outward = self._grad(location)
        else:
            sgn_dist, delta, *_ = self.approximate_closest_surface(location)
            outward = math.vec_normalize(math.sign(-sgn_dist) * delta)
        return sgn_dist, outward

    def approximate_signed_distance(self, location: Tensor) -> Tensor:
        return self._sdf(location)

    def sample_uniform(self, *shape: math.Shape) -> Tensor:
        raise NotImplementedError

    def bounding_radius(self) -> Tensor:
        return self._bounding_radius

    def bounding_half_extent(self) -> Tensor:
        return self._bounds.half_size  # this could be too small if the center is not in the middle of the bounds

    def bounding_box(self) -> 'BaseBox':
        return self._bounds

    def shifted(self, delta: Tensor) -> 'Geometry':
        raise NotImplementedError("SDF does not yet support shifting")

    def at(self, center: Tensor) -> 'Geometry':
        raise NotImplementedError("SDF does not yet support shifting")

    def rotated(self, angle: Union[float, Tensor]) -> 'Geometry':
        raise NotImplementedError("SDF does not yet support rotation")

    def scaled(self, factor: Union[float, Tensor]) -> 'Geometry':
        off_center = self._center - self._bounds.center
        volume = self._volume * factor ** self.spatial_rank
        bounds = self._bounds.scaled(factor).shifted(off_center * (factor - 1)).corner_representation()
        return SDF(self._sdf, bounds, self._center, volume, self._bounding_radius * factor)

    def __getitem__(self, item):
        item = slicing_dict(self, item)
        if not item:
            return self
        raise NotImplementedError

    @staticmethod
    def __stack__(values: tuple, dim: Shape, **kwargs) -> 'Geometry':
        from ._geom_ops import GeometryStack
        return GeometryStack(math.layout(values, dim))

Ancestors

phi.geom._geom.Geometry

Instance variables

prop boundary_elements : Dict[Any, Dict[str, slice]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[Any, Dict[str, slice]]:
    return {}

prop boundary_faces : Dict[Any, Dict[str, slice]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
    return {}

prop bounds : phi.geom._box.BaseBox

Expand source code

@property
def bounds(self) -> BaseBox:
    return self._bounds

prop center : phiml.math._tensors.Tensor

Center location in single channel dimension.

Expand source code

@property
def center(self) -> Tensor:
    return self._center

prop corners : phiml.math._tensors.Tensor

Returns

Expand source code

@property
def corners(self) -> Tensor:
    raise NotImplementedError(f"SDF does not support corners")

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    raise NotImplementedError(f"SDF does not support faces")

prop face_centers : phiml.math._tensors.Tensor

Expand source code

@property
def face_centers(self) -> Tensor:
    raise NotImplementedError(f"SDF does not support faces")

prop face_normals : phiml.math._tensors.Tensor

Normal vectors of cell faces, including boundary faces. Shape (elements, ~, vector). For meshes, The vectors point out of the primal cells and into the dual cells.

Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.

Expand source code

@property
def face_normals(self) -> Tensor:
    raise NotImplementedError(f"SDF does not support faces")

prop face_shape : phiml.math._shape.Shape

Returns

Expand source code

@property
def face_shape(self) -> Shape:
    return math.EMPTY_SHAPE

prop faces : Geometry

Expand source code

@property
def faces(self) -> 'Geometry':
    raise NotImplementedError(f"SDF does not support faces")

prop shape : phiml.math._shape.Shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self) -> Shape:
    return self._out_shape & self._bounds.shape

prop size

Expand source code

@property
def size(self):
    return self._bounds.size

prop volume : phiml.math._tensors.Tensor

phi.math.Tensor representing the volume of each element. The result retains batch, spatial and instance dimensions.

Expand source code

@property
def volume(self) -> Tensor:
    return self._volume

Methods

def approximate_closest_surface(self, location: phiml.math._tensors.Tensor, refine_iter=0) ‑> Tuple[phiml.math._tensors.Tensor, phiml.math._tensors.Tensor, phiml.math._tensors.Tensor, phiml.math._tensors.Tensor, phiml.math._tensors.Tensor]

Find the closest surface face of this geometry given a point that can be outside or inside the geometry.

Args

location: Tensor with a single channel dimension called vector. Can have arbitrary other dimensions.

Returns

signed_distance: Scalar signed distance from location to the closest point on the surface. Positive values indicate the point lies outside the geometry, negative values indicate the point lies inside the geometry.
delta: Vector-valued distance vector from location to the closest point on the surface.
normal: Closest surface normal vector.
offset: Min distance of a surface-tangential plane from 0 as a scalar.
face_index: (Optional) An index vector pointing at the closest face.

def approximate_signed_distance(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Computes the approximate distance from location to the surface of the geometry. Locations outside return positive values, inside negative values and zero exactly at the boundary.

When dealing with unions or collections of geometries (instance dimensions), the shortest distance to any instance is returned. This also holds for negative distances.

Args

location: Tensor with one channel dim vector matching the geometry's vector dim.

Returns

Float Tensor

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_box(self) ‑> phi.geom._box.BaseBox

Returns the approximately smallest axis-aligned box that contains this Geometry. The center of the box may not be equal to self.center.

Returns

Box or Cuboid that fully contains this Geometry.

def bounding_half_extent(self) ‑> phiml.math._tensors.Tensor

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def bounding_radius(self) ‑> phiml.math._tensors.Tensor

Returns the radius of a Sphere object that fully encloses this geometry. The sphere is centered at the center of this geometry.

def lies_inside(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def rotated(self, angle: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def sample_uniform(self, *shape: phiml.math._shape.Shape) ‑> phiml.math._tensors.Tensor

Samples uniformly distributed random points inside this volume.

Args

*shape: How many points to sample per individual geometry.

Returns

Tensor containing all dimensions from Geometry.shape, shape as well as a channel dimension vector matching the dimensionality of this Geometry.

def scaled(self, factor: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Scales each individual geometry by factor. The individual center points act as pivots for the operation.

Args

factor: Returns:

def sdf_and_gradient(self, location: phiml.math._tensors.Tensor, refine_iter=0) ‑> Tuple[phiml.math._tensors.Tensor, phiml.math._tensors.Tensor]

def shifted(self, delta: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a translated version of this geometry.

See Also: Geometry.at().

Args

delta: direction vector
delta: Tensor:

Returns

Geometry: shifted geometry

class SDFGrid (sdf: phiml.math._tensors.Tensor, bounds: phi.geom._box.BaseBox, approximate_outside=True, gradient: phiml.math._tensors.Tensor = None, to_surface: phiml.math._tensors.Tensor = None, surf_normal: phiml.math._tensors.Tensor = None, surf_index: phiml.math._tensors.Tensor = None, center: phiml.math._tensors.Tensor = None, volume: phiml.math._tensors.Tensor = None, bounding_radius: phiml.math._tensors.Tensor = None)

Grid-based signed distance field.

Args

sdf: Signed distance values. Tensor with spatial dimensions corresponding to the physical space. Each value samples the SDF value at the center of a virtual cell.
bounds: Grid limits. The bounds fully enclose all virtual cells.
approximate_outside: Whether queries outside the SDF grid should return approximate values. This requires additional computations.
gradient: (Optional) Pre-computed gradient grid. Will be computed otherwise.
center: (Optional) Geometry center point. Will be computed otherwise.
volume: (Optional) Geometry volume. Will be computed otherwise.
bounding_radius: (Optional) Geometry bounding radius around center. Will be computed otherwise.

Expand source code

class SDFGrid(Geometry):
    """
    Grid-based signed distance field.
    """
    def __init__(self,
                 sdf: Tensor,
                 bounds: BaseBox,
                 approximate_outside=True,
                 gradient: Tensor = None,
                 to_surface: Tensor = None, surf_normal: Tensor = None, surf_index: Tensor = None,
                 center: Tensor = None, volume: Tensor = None, bounding_radius: Tensor = None):
        """
        Args:
            sdf: Signed distance values. `Tensor` with spatial dimensions corresponding to the physical space.
                Each value samples the SDF value at the center of a virtual cell.
            bounds: Grid limits. The bounds fully enclose all virtual cells.
            approximate_outside: Whether queries outside the SDF grid should return approximate values. This requires additional computations.
            gradient: (Optional) Pre-computed gradient grid. Will be computed otherwise.
            center: (Optional) Geometry center point. Will be computed otherwise.
            volume: (Optional) Geometry volume. Will be computed otherwise.
            bounding_radius: (Optional) Geometry bounding radius around center. Will be computed otherwise.
        """
        super().__init__()
        self._sdf = sdf
        self._bounds = bounds
        self._approximate_outside = approximate_outside
        dx = bounds.size / spatial(sdf)
        if gradient is True:
            grad = math.spatial_gradient(sdf, dx=dx, difference='forward', padding=math.extrapolation.ZERO_GRADIENT, stack_dim=channel('vector'))
            self._grad = grad[{dim: slice(0, -1) for dim in spatial(sdf).names}]
        else:
            self._grad = gradient
        self._to_surface = to_surface
        self._surf_normal = surf_normal
        self._surf_index = surf_index
        if center is not None:
            self._center = center
        else:
            min_index = math.argmin(self._sdf, spatial, index_dim=channel('vector'))
            self._center = bounds.local_to_global(min_index / spatial(sdf))
        if volume is not None:
            self._volume = volume
        else:
            filled = math.sum(sdf < 0)
            self._volume = filled * math.prod(dx)
        if bounding_radius is not None:
            self._bounding_radius = bounding_radius
        else:
            points = UniformGrid(spatial(sdf), self._bounds).center
            dist = math.vec_length(points - self._center)
            dist = math.where(self._sdf <= 0, dist, 0)
            self._bounding_radius = math.max(dist)

    def __variable_attrs__(self):
        return '_sdf', '_bounds', '_grad', '_to_surface', '_surf_normal', '_surf_index', '_center', '_volume', '_bounding_radius'

    def __value_attrs__(self):
        return '_sdf',

    @property
    def values(self):
        """Signed distance grid."""
        return self._sdf

    def with_values(self, values: Tensor):
        values = expand(values, spatial(self._sdf) - spatial(values))
        return SDFGrid(values, self._bounds, self._approximate_outside, self._grad, self._center, self._volume, self._bounding_radius)

    @property
    def bounds(self) -> BaseBox:
        return self._bounds

    @property
    def size(self) -> Tensor:
        return self._bounds.size

    @property
    def resolution(self) -> Shape:
        return spatial(self._sdf)

    @property
    def dx(self) -> Tensor:
        return self._bounds.size / spatial(self._sdf)

    @property
    def points(self) -> Tensor:
        return UniformGrid(spatial(self._sdf), self._bounds).center

    @property
    def grid(self) -> UniformGrid:
        return UniformGrid(spatial(self._sdf), self._bounds)

    @property
    def center(self) -> Tensor:
        return self._center

    @property
    def shape(self) -> Shape:
        return non_spatial(self._sdf) & channel(vector=spatial(self._sdf))

    @property
    def volume(self) -> Tensor:
        return self._volume

    @property
    def faces(self) -> 'Geometry':
        raise NotImplementedError(f"SDF does not support faces")

    @property
    def face_centers(self) -> Tensor:
        raise NotImplementedError(f"SDF does not support faces")

    @property
    def face_areas(self) -> Tensor:
        raise NotImplementedError(f"SDF does not support faces")

    @property
    def face_normals(self) -> Tensor:
        raise NotImplementedError(f"SDF does not support faces")

    @property
    def boundary_elements(self) -> Dict[Any, Dict[str, slice]]:
        return {}

    @property
    def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
        return {}

    @property
    def face_shape(self) -> Shape:
        return math.EMPTY_SHAPE

    @property
    def corners(self) -> Tensor:
        raise NotImplementedError(f"SDF does not support corners")

    def lies_inside(self, location: Tensor) -> Tensor:
        float_idx = (location - self._bounds.lower) / self.size * self.resolution
        sdf_val = math.grid_sample(self._sdf, float_idx - .5, math.extrapolation.ZERO_GRADIENT)
        if self._approximate_outside:
            within_bounds = self._bounds.lies_inside(location)
            return within_bounds & (sdf_val <= 0)
        else:
            return sdf_val <= 0

    def approximate_closest_surface(self, location: Tensor) -> Tuple[Tensor, Tensor, Tensor, Tensor, Tensor]:
        if self._approximate_outside:
            location = self._bounds.push(location, outward=False)
        float_idx = (location - self._bounds.lower) / self.size * self.resolution
        sgn_dist = math.grid_sample(self._sdf, float_idx - .5, math.extrapolation.ZERO_GRADIENT)
        if self._to_surface is not None:
            to_surf = math.grid_sample(self._to_surface, float_idx - .5, math.extrapolation.ZERO_GRADIENT)
        else:
            sdf_grad = math.grid_sample(self._grad, float_idx - 1, math.extrapolation.ZERO_GRADIENT)
            sdf_grad = math.vec_normalize(sdf_grad)  # theoretically not necessary
            to_surf = sgn_dist * -sdf_grad
        surface_pos = location + to_surf
        if self._surf_normal is not None:
            normal = math.grid_sample(self._surf_normal, float_idx - .5, math.extrapolation.ZERO_GRADIENT)
            int_index = clip(math.to_int32(float_idx), 0, wrap(spatial(self._surf_index), '(x,y,z)')-1)
            face_index = self._surf_index[int_index]
        else:
            surf_float_idx = (surface_pos - self._bounds.lower) / self.size * self.resolution
            normal = math.grid_sample(self._grad, surf_float_idx - 1, math.extrapolation.ZERO_GRADIENT)
            # normal = math.where(self._bounds.lies_inside(surface_pos), normal, sdf_grad)  # use current normal if surface point is outside SDF grid
            normal = math.vec_normalize(normal)
            face_index = None
        offset = normal.vector @ surface_pos.vector
        return sgn_dist, to_surf, normal, offset, face_index

    def approximate_signed_distance(self, location: Tensor) -> Tensor:
        float_idx = (location - self._bounds.lower) / self.size * self.resolution
        sdf_val = math.grid_sample(self._sdf, float_idx - .5, math.extrapolation.ZERO_GRADIENT)
        if self._approximate_outside:
            within_bounds = self._bounds.lies_inside(location)
            dist_from_center = math.vec_length(location - self._center) - self._bounding_radius
            return math.where(within_bounds, sdf_val, dist_from_center)
        else:
            return sdf_val

    def sample_uniform(self, *shape: math.Shape) -> Tensor:
        raise NotImplementedError

    def bounding_radius(self) -> Tensor:
        return self._bounding_radius

    def bounding_half_extent(self) -> Tensor:
        return self._bounds.half_size  # this could be too small if the center is not in the middle of the bounds

    def shifted(self, delta: Tensor) -> 'Geometry':
        return SDFGrid(self._sdf, self._bounds.shifted(delta), self._approximate_outside, self._grad, self._to_surface, self._surf_normal, self._surf_index, self._center + delta, self._volume, self._bounding_radius)

    def at(self, center: Tensor) -> 'Geometry':
        return self.shifted(center - self._center)

    def rotated(self, angle: Union[float, Tensor]) -> 'Geometry':
        raise NotImplementedError("SDF does not yet support rotation")

    def scaled(self, factor: Union[float, Tensor]) -> 'Geometry':
        off_center = self._center - self._bounds.center
        volume = self._volume * factor ** self.spatial_rank
        bounds = self._bounds.scaled(factor).shifted(off_center * (factor - 1)).corner_representation()
        return SDFGrid(self._sdf, bounds, self._approximate_outside, self._grad, self._to_surface, self._surf_normal, self._surf_index, self._center, volume, self._bounding_radius * factor)

    def __getitem__(self, item):
        item = slicing_dict(self, item)
        if 'vector' in item:
            raise NotImplementedError("SDF projection not yet supported")
        return SDFGrid(self._sdf[item], self._bounds[item], self._approximate_outside, math.slice(self._grad, item), math.slice(self._to_surface, item), math.slice(self._surf_normal, item), math.slice(self._surf_index, item), math.slice(self._center, item), math.slice(self._volume, item), math.slice(self._bounding_radius, item))

    def approximate_occupancy(self):
        assert self._surf_normal is not None
        unit_corners = Cuboid(half_size=.5*self.dx).corners
        surf_dist = self._surf_normal.vector @ self._to_surface.vector
        corner_sdf = unit_corners.vector @ self._surf_normal.vector - surf_dist
        total_dist = math.sum(abs(corner_sdf), dual)
        neg_dist = math.sum(math.minimum(0, corner_sdf), dual)
        occ_near_surf = -neg_dist / total_dist
        occ_away = self._sdf <= 0
        return math.where(abs(self._sdf) < .5*math.vec_length(self.dx), occ_near_surf, occ_away)

    def approximate_fraction_inside(self, other_geometry: 'Geometry', balance: Union[Tensor, Number] = 0.5) -> Tensor:
        if other_geometry == self.grid and math.always_close(balance, .5):
            return self.approximate_occupancy()
        else:
            return Geometry.approximate_fraction_inside(self, other_geometry, balance)

    def downsample2x(self):
        s, g, t, n, i = [math.downsample2x(v) for v in (self._sdf, self._grad, self._to_surface, self._surf_normal, self._surf_index)]
        return SDFGrid(s, self._bounds, self._approximate_outside, g, t, n, i, self._center, self._volume, self._bounding_radius)

Ancestors

phi.geom._geom.Geometry

Instance variables

prop boundary_elements : Dict[Any, Dict[str, slice]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[Any, Dict[str, slice]]:
    return {}

prop boundary_faces : Dict[Any, Dict[str, slice]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
    return {}

prop bounds : phi.geom._box.BaseBox

Expand source code

@property
def bounds(self) -> BaseBox:
    return self._bounds

prop center : phiml.math._tensors.Tensor

Center location in single channel dimension.

Expand source code

@property
def center(self) -> Tensor:
    return self._center

prop corners : phiml.math._tensors.Tensor

Returns

Expand source code

@property
def corners(self) -> Tensor:
    raise NotImplementedError(f"SDF does not support corners")

prop dx : phiml.math._tensors.Tensor

Expand source code

@property
def dx(self) -> Tensor:
    return self._bounds.size / spatial(self._sdf)

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    raise NotImplementedError(f"SDF does not support faces")

prop face_centers : phiml.math._tensors.Tensor

Expand source code

@property
def face_centers(self) -> Tensor:
    raise NotImplementedError(f"SDF does not support faces")

prop face_normals : phiml.math._tensors.Tensor

Normal vectors of cell faces, including boundary faces. Shape (elements, ~, vector). For meshes, The vectors point out of the primal cells and into the dual cells.

Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.

Expand source code

@property
def face_normals(self) -> Tensor:
    raise NotImplementedError(f"SDF does not support faces")

prop face_shape : phiml.math._shape.Shape

Returns

Expand source code

@property
def face_shape(self) -> Shape:
    return math.EMPTY_SHAPE

prop faces : Geometry

Expand source code

@property
def faces(self) -> 'Geometry':
    raise NotImplementedError(f"SDF does not support faces")

prop grid : phi.geom._grid.UniformGrid

Expand source code

@property
def grid(self) -> UniformGrid:
    return UniformGrid(spatial(self._sdf), self._bounds)

prop points : phiml.math._tensors.Tensor

Expand source code

@property
def points(self) -> Tensor:
    return UniformGrid(spatial(self._sdf), self._bounds).center

prop resolution : phiml.math._shape.Shape

Expand source code

@property
def resolution(self) -> Shape:
    return spatial(self._sdf)

prop shape : phiml.math._shape.Shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self) -> Shape:
    return non_spatial(self._sdf) & channel(vector=spatial(self._sdf))

prop size : phiml.math._tensors.Tensor

Expand source code

@property
def size(self) -> Tensor:
    return self._bounds.size

prop values

Signed distance grid.

Expand source code

@property
def values(self):
    """Signed distance grid."""
    return self._sdf

prop volume : phiml.math._tensors.Tensor

phi.math.Tensor representing the volume of each element. The result retains batch, spatial and instance dimensions.

Expand source code

@property
def volume(self) -> Tensor:
    return self._volume

Methods

Find the closest surface face of this geometry given a point that can be outside or inside the geometry.

Args

location: Tensor with a single channel dimension called vector. Can have arbitrary other dimensions.

Returns

signed_distance: Scalar signed distance from location to the closest point on the surface. Positive values indicate the point lies outside the geometry, negative values indicate the point lies inside the geometry.
delta: Vector-valued distance vector from location to the closest point on the surface.
normal: Closest surface normal vector.
offset: Min distance of a surface-tangential plane from 0 as a scalar.
face_index: (Optional) An index vector pointing at the closest face.

def approximate_fraction_inside(self, other_geometry: Geometry, balance: phiml.math._tensors.Tensor | numbers.Number = 0.5) ‑> phiml.math._tensors.Tensor

To call this method on batches of geometries of same shape, pass a batched Geometry instance. The result tensor will match the batch shape of other_geometry.

The result may only be accurate in special cases. The given geometries may be approximated as spheres or boxes using bounding_radius() and bounding_half_extent().

The default implementation of this method approximates other_geometry as a Sphere and computes the fraction using approximate_signed_distance().

Args

other_geometry: Geometry or geometry batch for which to compute the overlap with self.
balance: Mid-level between 0 and 1, default 0.5. This value is returned when exactly half of other_geometry lies inside self. 0.5 < balance <= 1 makes self seem larger while 0 <= balance < 0.5makes self seem smaller.

Returns

fraction of cell volume lying inside the geometry. float tensor of shape (other_geometry.batch_shape, 1).

def approximate_occupancy(self)

def approximate_signed_distance(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Computes the approximate distance from location to the surface of the geometry. Locations outside return positive values, inside negative values and zero exactly at the boundary.

When dealing with unions or collections of geometries (instance dimensions), the shortest distance to any instance is returned. This also holds for negative distances.

Args

location: Tensor with one channel dim vector matching the geometry's vector dim.

Returns

Float Tensor

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_half_extent(self) ‑> phiml.math._tensors.Tensor

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def bounding_radius(self) ‑> phiml.math._tensors.Tensor

Returns the radius of a Sphere object that fully encloses this geometry. The sphere is centered at the center of this geometry.

def downsample2x(self)

def lies_inside(self, location: phiml.math._tensors.Tensor) ‑> phiml.math._tensors.Tensor

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def rotated(self, angle: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def sample_uniform(self, *shape: phiml.math._shape.Shape) ‑> phiml.math._tensors.Tensor

Samples uniformly distributed random points inside this volume.

Args

*shape: How many points to sample per individual geometry.

Returns

Tensor containing all dimensions from Geometry.shape, shape as well as a channel dimension vector matching the dimensionality of this Geometry.

def scaled(self, factor: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Scales each individual geometry by factor. The individual center points act as pivots for the operation.

Args

factor: Returns:

def shifted(self, delta: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a translated version of this geometry.

See Also: Geometry.at().

Args

delta: direction vector
delta: Tensor:

Returns

Geometry: shifted geometry

def with_values(self, values: phiml.math._tensors.Tensor)

class Sphere (center: phiml.math._tensors.Tensor = None, radius: phiml.math._tensors.Tensor | float = None, volume: phiml.math._tensors.Tensor | float = None, radius_variable=True, **center_: phiml.math._tensors.Tensor | float)

N-dimensional sphere. Defined through center position and radius.

Args

center: Sphere center as Tensor with vector dimension. The spatial dimension order should be specified in the vector dimension via item names. Can be left empty to specify dimensions via kwargs.
radius: Sphere radius as float or Tensor
**center_: Specifies center when the center argument is not given. Center position by dimension, e.g. x=0.5, y=0.2.

Expand source code

class Sphere(Geometry):
    """
    N-dimensional sphere.
    Defined through center position and radius.
    """

    def __init__(self,
                 center: Tensor = None,
                 radius: Union[float, Tensor] = None,
                 volume: Union[float, Tensor] = None,
                 radius_variable=True,
                 **center_: Union[float, Tensor]):
        """
        Args:
            center: Sphere center as `Tensor` with `vector` dimension.
                The spatial dimension order should be specified in the `vector` dimension via item names.
                Can be left empty to specify dimensions via kwargs.
            radius: Sphere radius as `float` or `Tensor`
            **center_: Specifies center when the `center` argument is not given. Center position by dimension, e.g. `x=0.5, y=0.2`.
        """
        if center is not None:
            assert isinstance(center, Tensor), f"center must be a Tensor but got {type(center).__name__}"
            assert 'vector' in center.shape, f"Sphere center must have a 'vector' dimension."
            assert center.shape.get_item_names('vector') is not None, f"Vector dimension must list spatial dimensions as item names. Use the syntax Sphere(x=x, y=y) to assign names."
            self._center = center
        else:
            self._center = wrap(tuple(center_.values()), math.channel(vector=tuple(center_.keys())))
        if radius is None:
            assert volume is not None, f"Either radius or volume must be specified but got neither."
            self._radius = Sphere.radius_from_volume(wrap(volume), self._center.vector.size)
        else:
            self._radius = wrap(radius)
        self._radius_variable = radius_variable
        assert 'vector' not in self._radius.shape, f"Sphere radius must not vary along vector but got {radius}"

    def __all_attrs__(self) -> tuple:
        return ('_center', '_radius')

    def __variable_attrs__(self) -> tuple:
        return ('_center', '_radius') if self._radius_variable else ('_center',)

    def __value_attrs__(self) -> tuple:
        return ()

    @property
    def shape(self):
        if self._center is None or self._radius is None:
            return None
        return self._center.shape & self._radius.shape

    @property
    def radius(self):
        return self._radius

    @property
    def center(self):
        return self._center

    @property
    def volume(self) -> math.Tensor:
        return Sphere.volume_from_radius(self._radius, self.spatial_rank)

    @staticmethod
    def volume_from_radius(radius: Union[float, Tensor], spatial_rank: int):
        if spatial_rank == 1:
            return 2 * radius
        elif spatial_rank == 2:
            return PI * radius ** 2
        elif spatial_rank == 3:
            return 4/3 * PI * radius ** 3
        else:
            raise NotImplementedError(f"spatial_rank>3 not supported, got {spatial_rank}")
            # n = self.spatial_rank
            # return math.pi ** (n // 2) / math.faculty(math.ceil(n / 2)) * self._radius ** n

    @staticmethod
    def radius_from_volume(volume: Union[float, Tensor], spatial_rank: int):
        if spatial_rank == 1:
            return volume / 2
        elif spatial_rank == 2:
            return math.sqrt(volume / PI)
        elif spatial_rank == 3:
            return (.75 / PI * volume) ** (1/3)
        else:
            raise NotImplementedError(f"spatial_rank>3 not supported, got {spatial_rank}")

    def lies_inside(self, location):
        distance_squared = math.sum((location - self.center) ** 2, dim='vector')
        return math.any(distance_squared <= self.radius ** 2, self.shape.instance)  # union for instance dimensions

    def approximate_signed_distance(self, location: Union[Tensor, tuple]):
        """
        Computes the exact distance from location to the closest point on the sphere.
        Very close to the sphere center, the distance takes a constant value.

        Args:
          location: float tensor of shape (batch_size, ..., rank)

        Returns:
          float tensor of shape (*location.shape[:-1], 1).

        """
        distance = math.vec_length(location - self._center, eps=1e-3)
        return math.min(distance - self.radius, self.shape.instance)  # union for instance dimensions

    def approximate_closest_surface(self, location: Tensor) -> Tuple[Tensor, Tensor, Tensor, Tensor, Tensor]:
        self_center = self.center
        self_radius = self.radius
        center_delta = location - self_center
        center_dist = math.vec_length(center_delta)
        sgn_dist = center_dist - self_radius
        if instance(self):
            self_center, self_radius, sgn_dist, center_delta, center_dist = math.at_min((self.center, self.radius, sgn_dist, center_delta, center_dist), key=abs(sgn_dist), dim=instance(self))
        normal = math.safe_div(center_delta, center_dist)
        default_normal = wrap([1] + [0] * (self.spatial_rank-1), self.shape['vector'])
        normal = math.where(center_dist == 0, default_normal, normal)
        surface_pos = self_center + self_radius * normal
        delta = surface_pos - location
        face_index = expand(0, non_channel(location))
        offset = normal.vector @ surface_pos.vector
        return sgn_dist, delta, normal, offset, face_index

    def sample_uniform(self, *shape: math.Shape):
        # --- Choose a distance from the center of the sphere, equally weighted by mass ---
        uniform = math.random_uniform(self.shape.non_singleton.without('vector'), *shape)
        if self.spatial_rank == 1:
            r = self.radius * uniform
        else:
            r = self.radius * (uniform ** (1 / self.spatial_rank))
        # --- Uniformly sample a unit vector for direction over the surface of the sphere (Muller 1959, Marsaglia 1972) ---
        unit_vector = math.random_normal(self.shape.non_singleton.without('vector'), *shape, self.shape['vector'])
        unit_vector /= math.vec_length(unit_vector, vec_dim='vector')
        return self.center + r * unit_vector

    def bounding_radius(self):
        return self.radius

    def bounding_half_extent(self):
        return expand(self.radius, self._center.shape.only('vector'))

    def at(self, center: Tensor) -> 'Geometry':
        return Sphere(center, self._radius, radius_variable=self._radius_variable)

    def rotated(self, angle):
        return self

    def scaled(self, factor: Union[float, Tensor]) -> 'Geometry':
        return Sphere(self.center, self.radius * factor, radius_variable=self._radius_variable)

    def __getitem__(self, item):
        item = slicing_dict(self, item)
        return Sphere(self._center[_keep_vector(item)], self._radius[item], radius_variable=self._radius_variable)

    @property
    def faces(self) -> 'Geometry':
        raise NotImplementedError(f"Sphere.faces not implemented.")

    @property
    def face_centers(self) -> Tensor:
        return math.zeros(self.shape & dual(shell=0))

    @property
    def face_areas(self) -> Tensor:
        return math.zeros(self.face_shape)

    @property
    def face_normals(self) -> Tensor:
        return math.zeros(self.shape & dual(shell=0))

    @property
    def boundary_elements(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
        return {}

    @property
    def boundary_faces(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
        return {}

    @property
    def face_shape(self) -> Shape:
        return self.shape.without('vector') & dual(shell=0)

    @property
    def corners(self) -> Tensor:
        return math.zeros(self.shape & dual(corners=0))

Ancestors

phi.geom._geom.Geometry

Static methods

def radius_from_volume(volume: phiml.math._tensors.Tensor | float, spatial_rank: int)
def volume_from_radius(radius: phiml.math._tensors.Tensor | float, spatial_rank: int)

Instance variables

prop boundary_elements : Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
    return {}

prop boundary_faces : Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[str, Tuple[Dict[str, slice], Dict[str, slice]]]:
    return {}

prop center

Center location in single channel dimension.

Expand source code

@property
def center(self):
    return self._center

prop corners : phiml.math._tensors.Tensor

Returns

Expand source code

@property
def corners(self) -> Tensor:
    return math.zeros(self.shape & dual(corners=0))

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    return math.zeros(self.face_shape)

prop face_centers : phiml.math._tensors.Tensor

Expand source code

@property
def face_centers(self) -> Tensor:
    return math.zeros(self.shape & dual(shell=0))

prop face_normals : phiml.math._tensors.Tensor

Normal vectors of cell faces, including boundary faces. Shape (elements, ~, vector). For meshes, The vectors point out of the primal cells and into the dual cells.

Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.

Expand source code

@property
def face_normals(self) -> Tensor:
    return math.zeros(self.shape & dual(shell=0))

prop face_shape : phiml.math._shape.Shape

Returns

Expand source code

@property
def face_shape(self) -> Shape:
    return self.shape.without('vector') & dual(shell=0)

prop faces : Geometry

Expand source code

@property
def faces(self) -> 'Geometry':
    raise NotImplementedError(f"Sphere.faces not implemented.")

prop radius

Expand source code

@property
def radius(self):
    return self._radius

prop shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self):
    if self._center is None or self._radius is None:
        return None
    return self._center.shape & self._radius.shape

prop volume : phiml.math._tensors.Tensor

phi.math.Tensor representing the volume of each element. The result retains batch, spatial and instance dimensions.

Expand source code

@property
def volume(self) -> math.Tensor:
    return Sphere.volume_from_radius(self._radius, self.spatial_rank)

Methods

Find the closest surface face of this geometry given a point that can be outside or inside the geometry.

Args

location: Tensor with a single channel dimension called vector. Can have arbitrary other dimensions.

Returns

signed_distance: Scalar signed distance from location to the closest point on the surface. Positive values indicate the point lies outside the geometry, negative values indicate the point lies inside the geometry.
delta: Vector-valued distance vector from location to the closest point on the surface.
normal: Closest surface normal vector.
offset: Min distance of a surface-tangential plane from 0 as a scalar.
face_index: (Optional) An index vector pointing at the closest face.

def approximate_signed_distance(self, location: phiml.math._tensors.Tensor | tuple)

Computes the exact distance from location to the closest point on the sphere. Very close to the sphere center, the distance takes a constant value.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

float tensor of shape (*location.shape[:-1], 1).

def at(self, center: phiml.math._tensors.Tensor) ‑> phi.geom._geom.Geometry

Returns a copy of this Geometry with the center at center. This is equal to calling self @ center.

See Also: Geometry.shifted().

Args

center: New center as Tensor.

Returns

Geometry.

def bounding_half_extent(self)

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def bounding_radius(self)

Returns the radius of a Sphere object that fully encloses this geometry. The sphere is centered at the center of this geometry.

def lies_inside(self, location)

Tests whether the given location lies inside or outside of the geometry. Locations on the surface count as inside.

When dealing with unions or collections of geometries (instance dimensions), a point lies inside the geometry if it lies inside any instance.

Args

location: float tensor of shape (batch_size, …, rank)

Returns

bool tensor of shape (*location.shape[:-1], 1).

def rotated(self, angle)

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def sample_uniform(self, *shape: phiml.math._shape.Shape)

Samples uniformly distributed random points inside this volume.

Args

*shape: How many points to sample per individual geometry.

Returns

Tensor containing all dimensions from Geometry.shape, shape as well as a channel dimension vector matching the dimensionality of this Geometry.

def scaled(self, factor: phiml.math._tensors.Tensor | float) ‑> phi.geom._geom.Geometry

Scales each individual geometry by factor. The individual center points act as pivots for the operation.

Args

factor: Returns:

class UniformGrid (resolution: phiml.math._shape.Shape = None, bounds: phi.geom._box.BaseBox = None, **resolution_)

An instance of UniformGrid represents all cells of a regular grid as a batch of boxes.

Expand source code

class UniformGrid(BaseBox):
    """
    An instance of UniformGrid represents all cells of a regular grid as a batch of boxes.
    """

    def __init__(self, resolution: Shape = None, bounds: BaseBox = None, **resolution_):
        assert resolution is None or resolution.is_uniform, f"spatial dimensions must form a uniform grid but got {resolution}"
        resolution = (resolution or EMPTY_SHAPE).spatial & spatial(**resolution_)
        bounds = _get_bounds(bounds, resolution)
        assert set(bounds.vector.item_names) == set(resolution.names)
        self._resolution = resolution.only(bounds.vector.item_names, reorder=True)  # reorder only
        self._bounds = bounds
        self._shape = self._resolution & bounds.shape.non_spatial
        staggered_shapes = [self._shape.spatial.with_dim_size(dim, self._shape.get_size(dim) + 1) for dim in self.vector.item_names]
        self._face_shape = shape_stack(dual(vector=self.vector.item_names), *staggered_shapes)

    @property
    def resolution(self):
        return self._resolution

    @property
    def bounds(self):
        return self._bounds

    @property
    def spatial_rank(self) -> int:
        return self._resolution.spatial_rank

    @property
    def center(self):
        local_coords = math.meshgrid(**{dim.name: math.linspace(0.5 / dim.size, 1 - 0.5 / dim.size, dim) for dim in self.resolution})
        points = self.bounds.local_to_global(local_coords)
        return points

    def position_of(self, voxel_index: Tensor):
        voxel_index = rename_dims(voxel_index, channel, 'vector')
        return self._bounds.lower + (voxel_index+.5) / self.resolution * self._bounds.size

    def voxel_at(self, location: Tensor, clamp=True):
        float_idx = (location - self._bounds.lower) / self._bounds.size * self.resolution
        index = math.to_int32(float_idx)
        if clamp:
            index = math.clip(index, 0, wrap(self.resolution, channel('vector'))-1)
        return index

    @property
    def boundary_elements(self) -> Dict[Any, Dict[str, slice]]:
        return {}

    @property
    def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
        result = {}
        for dim in self.vector.item_names:
            result[dim+'-'] = {'~vector': dim, dim: slice(1)}
            result[dim+'+'] = {'~vector': dim, dim: slice(-1, None)}
        return result

    @property
    def face_centers(self) -> Tensor:
        centers = [self.stagger(dim, True, True).center for dim in self.vector.item_names]
        return stack(centers, dual(vector=self.vector.item_names))

    @property
    def faces(self) -> Geometry:
        slices = [self.stagger(d, True, True) for d in self.resolution.names]
        return stack(slices, dual(vector=self.vector.item_names))

    @property
    def face_normals(self) -> Tensor:
        normals = [vec(**{d: float(d == dim) for d in self.vector.item_names}) for dim in self.vector.item_names]
        return stack(normals, dual(vector=self.vector.item_names))

    @property
    def face_areas(self) -> Tensor:
        areas = [math.prod(self.dx.vector[[d for d in self.vector.item_names if d != dim]], 'vector') for dim in self.vector.item_names]
        return stack(areas, dual(vector=self.vector.item_names))

    @property
    def face_shape(self) -> Shape:
        return self._face_shape

    def interior(self) -> 'Geometry':
        raise GeometryException("Regular grid does not have an interior")

    @property
    def grid_size(self):
        return self._bounds.size

    @property
    def size(self):
        return self.bounds.size / math.wrap(self.resolution.sizes)

    @property
    def dx(self):
        return self.bounds.size / self.resolution

    @property
    def lower(self):
        return self.center - self.half_size

    @property
    def upper(self):
        return self.center + self.half_size

    @property
    def half_size(self):
        return self.bounds.size / self.resolution.sizes / 2

    @property
    def rotation_matrix(self) -> Optional[Tensor]:
        return None

    def with_scaled_resolution(self, scale: float):
        return UniformGrid(self._resolution * scale, self._bounds)

    def __getitem__(self, item):
        item = slicing_dict(self, item)
        resolution = self._resolution.after_gather(item)
        bounds = self._bounds[{d: s for d, s in item.items() if d != 'vector'}]
        if 'vector' in item:
            resolution = resolution.only(item['vector'], reorder=True)
            bounds = bounds.vector[item['vector']]
        bounds = bounds.vector[resolution.name_list]
        dx = self.size
        for dim, selection in item.items():
            if dim in resolution:
                if isinstance(selection, slice):
                    start = selection.start or 0
                    if start < 0:
                        start += self.resolution.get_size(dim)
                    stop = selection.stop or self.resolution.get_size(dim)
                    if stop < 0:
                        stop += self.resolution.get_size(dim)
                    assert selection.step is None or selection.step == 1
                else:  # int slices are not contained in resolution anymore
                    raise ValueError(f"Illegal selection: {item}")
                dim_mask = math.wrap(self.resolution.mask(dim))
                lower = bounds.lower + start * dim_mask * dx
                upper = bounds.upper + (stop - self.resolution.get_size(dim)) * dim_mask * dx
                bounds = Box(lower, upper)
        return UniformGrid(resolution, bounds)

    def __pack_dims__(self, dims: Tuple[str, ...], packed_dim: Shape, pos: Optional[int], **kwargs) -> 'Cuboid':
        return math.pack_dims(self.center_representation(size_variable=False), dims, packed_dim, pos, **kwargs)

    @staticmethod
    def __stack__(values: tuple, dim: Shape, **kwargs) -> 'Geometry':
        from ._geom_ops import GeometryStack
        set_op = kwargs.get('set_op')
        return GeometryStack(math.layout(values, dim), set_op)

    def __replace_dims__(self, dims: Tuple[str, ...], new_dims: Shape, **kwargs) -> 'UniformGrid':
        resolution = math.rename_dims(self._resolution, dims, new_dims).spatial
        bounds = math.rename_dims(self._bounds, dims, new_dims, **kwargs)[resolution.name_list]
        return UniformGrid(resolution, bounds)

    def list_cells(self, dim_name):
        center = math.pack_dims(self.center, self._shape.spatial.names, dim_name)
        return Cuboid(center, self.half_size, size_variable=False)

    def stagger(self, dim: str, lower: bool, upper: bool):
        dim_mask = np.array(self.resolution.mask(dim))
        unit = self.bounds.size / self.resolution * dim_mask
        bounds = Box(self.bounds.lower + unit * (-0.5 if lower else 0.5), self.bounds.upper + unit * (0.5 if upper else -0.5))
        ext_res = self.resolution.sizes + dim_mask * (int(lower) + int(upper) - 1)
        return UniformGrid(self.resolution.with_sizes(ext_res), bounds)

    def staggered_cells(self, boundaries: Extrapolation) -> Dict[str, 'UniformGrid']:
        grids = {}
        for dim in self.vector.item_names:
            grids[dim] = self.stagger(dim, *boundaries.valid_outer_faces(dim))
        return grids

    def padded(self, widths: dict):
        resolution, bounds = self.resolution, self.bounds
        for dim, (lower, upper) in widths.items():
            masked_dx = self.dx * math.dim_mask(self.resolution, dim)
            resolution = resolution.with_dim_size(dim, self.resolution.get_size(dim) + lower + upper)
            bounds = Box(bounds.lower - masked_dx * lower, bounds.upper + masked_dx * upper)
        return UniformGrid(resolution, bounds)

    @property
    def shape(self):
        return self._shape

    def shifted(self, delta: Tensor, **delta_by_dim) -> BaseBox:
        # delta += math.padded_stack()
        if delta.shape.spatial_rank == 0:
            return UniformGrid(self.resolution, self.bounds.shifted(delta))
        else:
            center = self.center + delta
            return Cuboid(center, self.half_size, size_variable=False)

    def rotated(self, angle) -> Geometry:
        raise NotImplementedError("Grids cannot be rotated. Use center_representation() to convert it to Cuboids first.")

    def shallow_equals(self, other):
        return self == other

    def __repr__(self):
        return f"{self._resolution}, bounds={self._bounds}"

    def __variable_attrs__(self):
        return ()

    def __value_attrs__(self):
        return ()

    def __eq__(self, other):
        if not isinstance(other, UniformGrid):
            return False
        return self._resolution == other._resolution and self._bounds == other._bounds

    def __hash__(self):
        return hash(self._resolution) + hash(self._bounds)

    @property
    def _center(self):
        return self.center

    @property
    def _half_size(self):
        return self.half_size

    @property
    def normal(self) -> Tensor:
        raise GeometryException("UniformGrid does not have normals")

    def bounding_half_extent(self) -> Tensor:
        return self.half_size

Ancestors

phi.geom._box.BaseBox
phi.geom._geom.Geometry

Instance variables

prop boundary_elements : Dict[Any, Dict[str, slice]]

Slices on the primal dimensions to mark boundary elements. Grids and meshes have no boundary elements and return {}. Dynamic graphs can define boundary elements for obstacles and walls.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_elements(self) -> Dict[Any, Dict[str, slice]]:
    return {}

prop boundary_faces : Dict[Any, Dict[str, slice]]

Slices on the dual dimensions to mark boundary faces.

Regular grids use the keys (dim, is_upper) to identify boundaries. Unstructured meshes use string identifiers for the boundaries. Dynamic graphs return slices along the dual dimensions.

Returns

Map from name to slicing dict.

Expand source code

@property
def boundary_faces(self) -> Dict[Any, Dict[str, slice]]:
    result = {}
    for dim in self.vector.item_names:
        result[dim+'-'] = {'~vector': dim, dim: slice(1)}
        result[dim+'+'] = {'~vector': dim, dim: slice(-1, None)}
    return result

prop bounds

Expand source code

@property
def bounds(self):
    return self._bounds

prop center

Center location in single channel dimension.

Expand source code

@property
def center(self):
    local_coords = math.meshgrid(**{dim.name: math.linspace(0.5 / dim.size, 1 - 0.5 / dim.size, dim) for dim in self.resolution})
    points = self.bounds.local_to_global(local_coords)
    return points

prop dx

Expand source code

@property
def dx(self):
    return self.bounds.size / self.resolution

prop face_areas : phiml.math._tensors.Tensor

Area of face connecting a pair of cells. Shape (elements, ~). Returns 0 for unconnected cells.

Expand source code

@property
def face_areas(self) -> Tensor:
    areas = [math.prod(self.dx.vector[[d for d in self.vector.item_names if d != dim]], 'vector') for dim in self.vector.item_names]
    return stack(areas, dual(vector=self.vector.item_names))

prop face_centers : phiml.math._tensors.Tensor

Expand source code

@property
def face_centers(self) -> Tensor:
    centers = [self.stagger(dim, True, True).center for dim in self.vector.item_names]
    return stack(centers, dual(vector=self.vector.item_names))

prop face_normals : phiml.math._tensors.Tensor

Normal vectors of cell faces, including boundary faces. Shape (elements, ~, vector). For meshes, The vectors point out of the primal cells and into the dual cells.

Instance/spatial dimensions along which the normal does not vary may not be included in the result tensor's shape.

Expand source code

@property
def face_normals(self) -> Tensor:
    normals = [vec(**{d: float(d == dim) for d in self.vector.item_names}) for dim in self.vector.item_names]
    return stack(normals, dual(vector=self.vector.item_names))

prop face_shape : phiml.math._shape.Shape

Returns

Expand source code

@property
def face_shape(self) -> Shape:
    return self._face_shape

prop faces : phi.geom._geom.Geometry

Expand source code

@property
def faces(self) -> Geometry:
    slices = [self.stagger(d, True, True) for d in self.resolution.names]
    return stack(slices, dual(vector=self.vector.item_names))

prop grid_size

Expand source code

@property
def grid_size(self):
    return self._bounds.size

prop half_size

Expand source code

@property
def half_size(self):
    return self.bounds.size / self.resolution.sizes / 2

prop lower

Expand source code

@property
def lower(self):
    return self.center - self.half_size

prop normal : phiml.math._tensors.Tensor

Expand source code

@property
def normal(self) -> Tensor:
    raise GeometryException("UniformGrid does not have normals")

prop resolution

Expand source code

@property
def resolution(self):
    return self._resolution

prop rotation_matrix : phiml.math._tensors.Tensor | None

Expand source code

@property
def rotation_matrix(self) -> Optional[Tensor]:
    return None

prop shape

The shape of a Geometry consists of the following dimensions:

A single channel dimension called 'vector' specifying the physical space
Instance dimensions denote that this geometry consists of multiple copies in the same space
Spatial dimensions denote a crystal (repeating structure) of this geometric primitive in space
Batch dimensions indicate non-interacting versions of this geometry for parallelization only.

Expand source code

@property
def shape(self):
    return self._shape

prop size

Expand source code

@property
def size(self):
    return self.bounds.size / math.wrap(self.resolution.sizes)

prop spatial_rank : int

Number of spatial dimensions of the geometry, 1 = 1D, 2 = 2D, 3 = 3D, etc.

Expand source code

@property
def spatial_rank(self) -> int:
    return self._resolution.spatial_rank

prop upper

Expand source code

@property
def upper(self):
    return self.center + self.half_size

Methods

def bounding_half_extent(self) ‑> phiml.math._tensors.Tensor

The bounding half-extent sets a limit on the outer-most point for each coordinate axis. Each component is non-negative.

def interior(self) ‑> phi.geom._geom.Geometry

def list_cells(self, dim_name)

def padded(self, widths: dict)

def position_of(self, voxel_index: phiml.math._tensors.Tensor)

def rotated(self, angle) ‑> phi.geom._geom.Geometry

Returns a rotated version of this geometry. The geometry is rotated about its center point.

Args

angle

Delta rotation. Either

Angle(s): scalar angle in 2d or euler angles along vector in 3D or higher.
Matrix: d⨯d rotation matrix

Returns

Rotated Geometry

def shallow_equals(self, other)

Quick equality check. May return False even if other == self. However, if True is returned, the geometries are guaranteed to be equal.

The shallow_equals() check does not compare all tensor elements but merely checks whether the same tensors are referenced.

def shifted(self, delta: phiml.math._tensors.Tensor, **delta_by_dim) ‑> phi.geom._box.BaseBox

Returns a translated version of this geometry.

See Also: Geometry.at().

Args

delta: direction vector
delta: Tensor:

Returns

Geometry: shifted geometry

def stagger(self, dim: str, lower: bool, upper: bool)

def staggered_cells(self, boundaries: phiml.math.extrapolation.Extrapolation) ‑> Dict[str, phi.geom._grid.UniformGrid]

def voxel_at(self, location: phiml.math._tensors.Tensor, clamp=True)

def with_scaled_resolution(self, scale: float)