Source code for Muscat.Containers.Filters.FilterBase

# -*- coding: utf-8 -*-
#
# This file is subject to the terms and conditions defined in
# file 'LICENSE.txt', which is part of this source code package.
#
from __future__ import annotations
from typing import Optional, List, Callable, Any, Union, NamedTuple
from functools import reduce
from collections.abc import Iterable
import numpy as np

from Muscat.Types import ArrayLike, MuscatIndex
from Muscat.Containers.Mesh import Mesh
from Muscat.Containers.ElementsContainers import ElementContainerLike

from Muscat.Containers.Tags import Tags
import Muscat.Containers.ElementsDescription as ED

try:
    from Muscat.Containers.Filters.NativeFilters import IntersectOrdered1D, UnionOrdered1D
except:
    IntersectOrdered1D = lambda x, y : np.intersect1d(x, y, assume_unique=True)
    UnionOrdered1D = np.union1d
    print("Warning! Using python version for intersect1d and union1D (slow)")

[docs]class ElementsSelection(NamedTuple): """Class to store a selection of elements of the same type. This class is immutable. The public attributes of this class - self.mesh: Mesh = The selection mesh - self.indices: np.ndarray = The indices of the selected elements locally to elements - self.elements: ElementsContainer = The current elements of mesh - self.meshOffset: int = Number of elements on the mesh before this ElementsContainer - self.selectionOffset: int = Number of elements on the selection before this selection - self.elementType: str = Equivalent of self.elements.elementType - self.part: int = Part number (non 0 for partitioned meshes) """ mesh:Mesh indices: np.ndarray elements: ElementContainerLike meshOffset:int selectionOffset:int part:int = 0
[docs] def Size(self) -> int: """Return the size of the selection Returns ------- int size of self.indices """ return self.indices.shape[0]
@property def elementType(self) -> ED.ElementType: """Return the element type for the current selection Returns ------- ED.ElementType _description_ """ return self.elements.elementType
[docs] def GetSelectionSlice(self) -> slice: """Return a slice on a vector of length self.size(). .. code-block:: python res = np.empty(elementToBeExtracted, dtype=MuscaFloat) for selection in filter(mesh): res[selection.GetSelectionSlice()] = mesh.elemFields["energy"][selection.GetMeshSlice()] Returns ------- slice Slice to be used on a vector of size : the number of elements extracted by the filter """ return slice(self.selectionOffset, self.Size()+self.selectionOffset)
[docs] def GetMeshSlice(self) -> np.ndarray: """Return a slice on a vector of length mesh.GetNumberOfElements() .. code-block:: python res = np.empty(elementToBeExtracted, dtype=MuscaFloat) for selection in filter(mesh): res[selection.GetSelectionSlice()] = mesh.elemFields["energy"][selection.GetMeshSlice()] Returns ------- slice, nd.ndarray Slice to be used on a vector of size : the number of elements on the mesh """ return self.meshOffset + self.indices
[docs] def GetNumpyMask(self) -> np.ndarray: """Return this selection as a numpy mask This means the array is of size this.elements.GetNumberOfElements() and is false for the selected elements Returns ------- np.ndarray the numpy mask (true for masked elements, e.g. not selected elements) """ mask = np.ones(self.elements.GetNumberOfElements(), dtype=bool) mask[self.indices] = False return mask
def __eq__(self, other: Union[ElementsSelection, Any]) -> bool: """_summary_ Parameters ---------- other : Union[ElementsSelection, Any] _description_ Returns ------- bool _description_ """ if isinstance(other, ElementsSelection): if self is other: return True else: if self.mesh is not other.mesh: return False if not np.array_equal(self.indices, other.indices): return False if self.elements is not other.elements: return False if self.meshOffset != other.meshOffset: return False if self.selectionOffset != other.selectionOffset: return False if self.part != other.part: return False else: return False return True def __ne__(self, other: Union[ElementsSelection, Any]): return not (self == other) def __iter__(self): """"for Backward compatibility with BasicTools filter name, data, ids = ElementsSelection() or for name, data, ids in [Element/Node]Filter(myMesh): ... """ yield self.elementType yield self.elements yield self.indices
ZoneType = Callable[[ArrayLike], np.ndarray]
[docs]class FilterBase(): """Base class to construct node and element filters Parameters ---------- zones : Optional[Union[List[ZoneType], ZoneType] ], optional a ZoneType or a list of ZoneType (ImplicitGeometry for example), by default None eMask : Optional[ArrayLike], optional a boolean vector of the size of the object to filter, by default None nMask : Optional[ArrayLike], optional a boolean vector of the size of the object to filter, by default None etags : Optional[Union[List[str], str] ], optional An element tag name or a list of element tags to be used, by default None nTags : Optional[Union[List[str], str]], optional An nodal tag name of a list of nodal tags names to use, by default None eFilter : Optional[Union[ElementFilter, FilterOperatorBase] ] And extra element filter, by default none nFilter : Optional[Union[NodeFilter, FilterOperatorBase] ] And extra nodal filter, by default none """ def __init__(self, zone: Optional[Union[ZoneType, List[ZoneType]]] = None, eMask: Optional[ArrayLike] = None, nMask: Optional[ArrayLike] = None, eTag: Optional[Union[str, List[str]]] = None, nTag: Optional[Union[str, List[str]]] = None, eFilter: Optional[FilterBase] = None, nFilter: Optional[FilterBase] = None ) -> None: super().__init__() ### zones ### self.zones: List = list() if zone is not None: self.SetZones(zone) ### eMask ## self.eMask = None if eMask is not None: eMask = np.asarray(eMask, dtype=bool) if len(eMask.shape) > 1 or len(eMask.shape) < 1: raise RuntimeError("eMask must have exactly one dimension") self.eMask = eMask ### nMask ## self.nMask: Union[None,np.ndarray] = None if nMask is not None: nMask = np.asarray(nMask, dtype=bool) if len(nMask.shape) > 1 or len(nMask.shape) < 1: raise RuntimeError("nMask must have exactly one dimension") self.nMask = nMask ### eTags ### self.eTags: List = list() if eTag is not None: self.SetETags(eTag) ### nTags ### self.nTags: List = list() if nTag is not None: self.SetNTags(nTag) self.eFilter = eFilter self.nFilter = nFilter
[docs] def IsEquivalent(self, other: Any) -> bool: """To check if 2 element filter are equivalent Parameters ---------- other : Any other object to check the equivalency Returns ------- bool True if the two filters are equal """ if id(self) == id(other): return True if isinstance(other, self.__class__): if sorted(self.eTags) != sorted(other.eTags): return False if sorted(self.nTags) != sorted(other.nTags): return False if self.zones != other.zones: return False if not np.array_equal(self.eMask, other.eMask): return False if not np.array_equal(self.nMask, other.nMask): return False if self.eFilter != other.eFilter: return False if self.nFilter != other.nFilter: return False return True else: return False
[docs] def SetETags(self, tagNames: Union[List[str], str]) -> None: """Set the element tag to use this will erase the previous list Parameters ---------- tagNames : List[str] the list of string with the elementary tag names """ if isinstance(tagNames, str): return self.SetETags([tagNames]) self.eTags = list(tagNames)
[docs] def SetNTags(self, tagNames: Union[List[str], str]) -> None: """Set the nodal tag to use this will erase the previous list Parameters ---------- tagNames : List[str] the list of string with the nodal tag names """ if isinstance(tagNames, str): return self.SetNTags([tagNames]) self.nTags = list(tagNames)
[docs] def AddETag(self, tagName: str) -> None: """Add a tagname to the list of elementary tags to treat Parameters ---------- tagName : str the tag name to be added """ self.eTags.append(tagName)
[docs] def AddNTag(self, tagName: str) -> None: """Add a tagname to the list of nodal tags to treat Parameters ---------- tagName : str the tag name to be added """ self.nTags.append(tagName)
[docs] def SetZones(self, zonesList: Union[ZoneType, List[ZoneType]]) -> None: """Set the zone list to use this will erase the previous list Parameters ---------- zonesList : List[ZoneType] The list of zones to be used """ if not isinstance(zonesList, Iterable): return self.SetZones([zonesList]) self.zones = zonesList
[docs] def AddZone(self, zone: ZoneType) -> None: """Add a zone to the list of zones to be treated by the filter Parameters ---------- zone : ZoneType a callable object capable of taking one argument with the points positions, and returning a vector of size pos.shape[0] with negative values for the entities to be selected by the filter """ self.zones.append(zone)
def _CheckTags_(self, tags: Tags, tagNames: List[str]) -> Union[np.ndarray, None]: """Internal function to compute the indices to be treated based on the tags Parameters ---------- tags : _type_ The tags container numberOfObjects : _type_ the total number of object (number of points or number of element in the current element container) Returns ------- np.ndarray or None np.ndarray of the the indices to be kept None if no tags present on to filter """ if len(tagNames) == 0: return None tagList = list(tags[t] for t in tagNames if t in tags.keys()) if len(tagList) == 0: return np.empty((0), dtype=MuscatIndex) # fast path if len(tagList) == 1: return tagList[0].GetIds() return reduce(UnionOrdered1D, (tag.GetIds() for tag in tagList)) def _CheckMask_(self, mask: Union[np.ndarray, None], start: MuscatIndex, size: MuscatIndex) -> Union[np.ndarray, None]: """Internal function to compute the indices based on the mask Parameters ---------- mask : np.ndarray the mask to work on start : MuscatIndex the start position for the current request size : MuscatIndex the stop position for the current request Returns ------- np.ndarray or None np.ndarray of the the indices to be kept None if no mask present on to filter """ if mask is not None: return np.where(mask[start:start+size])[0].astype(MuscatIndex) return None def _Done_(self, data: Union[np.ndarray, None]) -> bool: """Helper function to check if we have to continue with the other tests Parameters ---------- data : Union[np.ndarray|None] the to check Returns ------- bool True if Data is not None and empty ( len == 0) False otherwise """ if data is not None and len(data) == 0: return True return False
[docs] def GetNodesIndices(self, mesh: Mesh) -> np.ndarray: """Get the nodes indices selected by this filter Parameters ---------- mesh : Mesh mesh to apply the filter Returns ------- Union np.ndarray The filtered nodes indices """ raise NotImplementedError("Cant call GetNodesIndices on FilterBase")
[docs]def Intersect1D(first: Union[np.ndarray, None], second: Union[np.ndarray, None]) -> Union[np.ndarray, None]: """Function to generate an intersection of two vectors (like np.intersect1d) but with the particularity of treat the case where the inputs can be None None represent a non filtered inputs, Parameters ---------- first : np.ndarray | None first vector of indices second : np.ndarray | None second vector of indices Returns ------- Union[np.ndarray,None] The intersection of two list """ if first is None: if second is None: return None else: return second else: if second is None: return first else: return IntersectOrdered1D(first,second)
[docs]def CheckIntegrity(GUI: bool = False) -> str: from Muscat.Helpers.CheckTools import MustFailFunction, MustFailFunctionWith from Muscat.Containers.Tags import Tags from Muscat.Containers.Filters.FilterObjects import ElementFilter, NodeFilter from Muscat.Containers import MeshCreationTools as UMCT import Muscat.Containers.ElementsDescription as ED mesh = UMCT.CreateSquare(dimensions=[10, 10], origin=[0, 0], spacing=[10/9, 10/9], ofTriangles=True) mesh2 = UMCT.CreateSquare(dimensions=[10, 10], origin=[0, 0], spacing=[10/9, 10/9], ofTriangles=True) es = ElementsSelection(mesh=mesh, indices=np.arange(5), elements=mesh.GetElementsOfType(ED.Bar_2), meshOffset=2, selectionOffset=3) np.testing.assert_equal(es.Size(), 5) np.testing.assert_equal(es.GetSelectionSlice(), slice(3, 8)) np.testing.assert_equal(es.GetMeshSlice(), [2, 3, 4, 5, 6]) es.GetNumpyMask() assert (es == es) es2 = ElementsSelection(mesh=mesh, indices=np.arange(5), elements=mesh.GetElementsOfType(ED.Bar_2), meshOffset=2, selectionOffset=3) assert (es == es2) es2 = ElementsSelection(mesh=mesh2, indices=np.arange(5), elements=mesh.GetElementsOfType(ED.Bar_2), meshOffset=2, selectionOffset=3) assert (es != es2) es2 = ElementsSelection(mesh=mesh, indices=np.arange(6), elements=mesh.GetElementsOfType(ED.Bar_2), meshOffset=2, selectionOffset=3) assert (es != es2) es2 = ElementsSelection(mesh=mesh, indices=np.arange(5), elements=mesh.GetElementsOfType(ED.Triangle_3), meshOffset=2, selectionOffset=3) assert (es != es2) es2 = ElementsSelection(mesh=mesh, indices=np.arange(5), elements=mesh.GetElementsOfType(ED.Bar_2), meshOffset=3, selectionOffset=3) assert (es != es2) es2 = ElementsSelection(mesh=mesh, indices=np.arange(5), elements=mesh.GetElementsOfType(ED.Bar_2), meshOffset=2, selectionOffset=4) assert (es != es2) es2 = ElementsSelection(mesh=mesh, indices=np.arange(5), elements=mesh.GetElementsOfType(ED.Bar_2), meshOffset=2, selectionOffset=3, part=1) assert (es != es2) assert (es != "some other type") a = np.array([1, 2]) b = np.array([2, 3]) np.testing.assert_equal(Intersect1D(None, b), b) np.testing.assert_equal(Intersect1D(a, None), a) np.testing.assert_equal(Intersect1D(a, b), [2]) np.testing.assert_equal(Intersect1D(None, None), None) obj = FilterBase(zone=[lambda x:x[:, 0] < 0, lambda x:x[:, 1] < 0], eMask=np.ones(2, dtype=MuscatIndex), nMask=np.ones(3, dtype=MuscatIndex), eTag=["A", "B"], nTag=["D", "E"], eFilter=ElementFilter(), nFilter=NodeFilter() ) obj.AddETag("C") obj.AddNTag("F") obj.AddZone(lambda x: x[:, 2] < 0) np.testing.assert_equal(obj._Done_(None), False) np.testing.assert_equal(obj._Done_(np.arange(0)), True) assert (obj.IsEquivalent(obj)) assert (obj.IsEquivalent(0) == False) obj2 = FilterBase() assert (obj.IsEquivalent(obj2) == False) obj2.SetETags(obj.eTags) assert (obj.IsEquivalent(obj2) == False) obj2.SetNTags(obj.nTags) assert (obj.IsEquivalent(obj2) == False) obj2.SetZones(obj.zones) assert (obj.IsEquivalent(obj2) == False) obj2.eMask = obj.eMask assert (obj.IsEquivalent(obj2) == False) obj2.nMask = obj.nMask assert (obj.IsEquivalent(obj2) == False) obj2.eFilter = obj.eFilter assert (obj.IsEquivalent(obj2) == False) obj2.nFilter = obj.nFilter assert (obj.IsEquivalent(obj2)) np.testing.assert_equal(obj._CheckMask_(obj.eMask, start=0, size=2), [0, 1]) np.testing.assert_equal(obj._CheckMask_(obj.nMask, start=0, size=3), [0, 1, 2]) np.testing.assert_equal(obj._CheckMask_(None, start=0, size=3), None) tags = Tags() tags.CreateTag("A").SetIds([0, 1]) tags.CreateTag("B").SetIds([1]) tags.CreateTag("C") np.testing.assert_equal(obj._CheckTags_(tags, obj.eTags), [0, 1]) np.testing.assert_equal(obj._CheckTags_(tags, ["B"]), [1]) np.testing.assert_equal(obj._CheckTags_(tags, []), None) np.testing.assert_equal(obj._CheckTags_(tags, ["X"]), []) FilterBase(eTag="Hello") FilterBase(nTag="Hello") FilterBase(zone="Hello") MustFailFunction(FilterBase, eMask=np.empty((2, 3))) MustFailFunction(FilterBase, nMask=np.empty((2, 3))) MustFailFunctionWith(NotImplementedError, FilterBase().GetNodesIndices,mesh) return "ok"
if __name__ == '__main__': print(CheckIntegrity(GUI=True)) # pragma: no cover