From 6937127e4772dc5f3fd16791ea343ecac764a191 Mon Sep 17 00:00:00 2001 From: Connor Ward Date: Fri, 31 Jan 2025 13:28:31 +0000 Subject: [PATCH 1/5] Make hashing of quadrature rules safe Also remove some uses of deprecated `@abstractproperty` and gives cells and point sets `repr()`s. --- FIAT/reference_element.py | 8 ++++++- finat/point_set.py | 32 ++++++++++++++++++++++---- finat/quadrature.py | 25 +++++++++++++++++--- test/finat/test_create_fiat_element.py | 5 ++++ test/finat/test_quadrature.py | 19 +++++++++++++++ 5 files changed, 80 insertions(+), 9 deletions(-) create mode 100644 test/finat/test_quadrature.py diff --git a/FIAT/reference_element.py b/FIAT/reference_element.py index 30fa94bb..af5c3bf1 100644 --- a/FIAT/reference_element.py +++ b/FIAT/reference_element.py @@ -126,7 +126,7 @@ def linalg_subspace_intersection(A, B): return U[:, :rank_c] -class Cell(object): +class Cell: """Abstract class for a reference cell. Provides accessors for geometry (vertex coordinates) as well as topology (orderings of vertices that make up edges, faces, etc.""" @@ -184,6 +184,9 @@ def __init__(self, shape, vertices, topology): # Dictionary with derived cells self._split_cache = {} + def __repr__(self): + return f"{type(self).__name__}({self.shape!r}, {self.vertices!r}, {self.topology!r})" + def _key(self): """Hashable object key data (excluding type).""" # Default: only type matters @@ -1130,6 +1133,9 @@ def __init__(self, *cells): super().__init__(TENSORPRODUCT, vertices, topology) self.cells = tuple(cells) + def __repr__(self): + return f"{type(self).__name__}({self.cells!r})" + def _key(self): return self.cells diff --git a/finat/point_set.py b/finat/point_set.py index 1497308c..5bb7f69e 100644 --- a/finat/point_set.py +++ b/finat/point_set.py @@ -1,4 +1,5 @@ -from abc import ABCMeta, abstractproperty +import abc +import hashlib from itertools import chain, product import numpy @@ -7,7 +8,7 @@ from gem.utils import cached_property -class AbstractPointSet(metaclass=ABCMeta): +class AbstractPointSet(abc.ABC): """A way of specifying a known set of points, perhaps with some (tensor) structure. @@ -15,8 +16,15 @@ class AbstractPointSet(metaclass=ABCMeta): where point_set_shape is () for scalar, (N,) for N element vector, (N, M) for N x M matrix etc. """ + def __hash__(self): + return int.from_bytes(hashlib.md5(repr(self).encode()).digest(), byteorder="big") - @abstractproperty + @abc.abstractmethod + def __repr__(self): + pass + + @property + @abc.abstractmethod def points(self): """A flattened numpy array of points or ``UnknownPointsArray`` object with shape (# of points, point dimension).""" @@ -27,12 +35,14 @@ def dimension(self): _, dim = self.points.shape return dim - @abstractproperty + @property + @abc.abstractmethod def indices(self): """GEM indices with matching shape and extent to the structure of the point set.""" - @abstractproperty + @property + @abc.abstractmethod def expression(self): """GEM expression describing the points, with free indices ``self.indices`` and shape (point dimension,).""" @@ -53,6 +63,9 @@ def __init__(self, point): assert len(point.shape) == 1 self.point = point + def __repr__(self): + return f"{type(self).__name__}({self.point!r})" + @cached_property def points(self): # Make sure we conform to the expected (# of points, point dimension) @@ -106,6 +119,9 @@ def __init__(self, points_expr): assert len(points_expr.shape) == 2 self._points_expr = points_expr + def __repr__(self): + return f"{type(self).__name__}({self._points_expr!r})" + @cached_property def points(self): return UnknownPointsArray(self._points_expr.shape) @@ -133,6 +149,9 @@ def __init__(self, points): assert len(points.shape) == 2 self.points = points + def __repr__(self): + return f"{type(self).__name__}({self.points!r})" + @cached_property def points(self): pass # set at initialisation @@ -177,6 +196,9 @@ class TensorPointSet(AbstractPointSet): def __init__(self, factors): self.factors = tuple(factors) + def __repr__(self): + return f"{type(self).__name__}({self.factors!r})" + @cached_property def points(self): return numpy.array([list(chain(*pt_tuple)) diff --git a/finat/quadrature.py b/finat/quadrature.py index ec6def12..6ea800a1 100644 --- a/finat/quadrature.py +++ b/finat/quadrature.py @@ -1,4 +1,5 @@ -from abc import ABCMeta, abstractproperty +import hashlib +from abc import ABCMeta, abstractmethod from functools import reduce import gem @@ -60,11 +61,23 @@ class AbstractQuadratureRule(metaclass=ABCMeta): """Abstract class representing a quadrature rule as point set and a corresponding set of weights.""" - @abstractproperty + def __hash__(self): + return int.from_bytes(hashlib.md5(repr(self).encode()).digest(), byteorder="big") + + def __eq__(self, other): + return type(other) is type(self) and hash(other) == hash(self) + + @abstractmethod + def __repr__(self): + pass + + @property + @abstractmethod def point_set(self): """Point set object representing the quadrature points.""" - @abstractproperty + @property + @abstractmethod def weight_expression(self): """GEM expression describing the weights, with the same free indices as the point set.""" @@ -110,6 +123,9 @@ def __init__(self, point_set, weights, ref_el=None, io_ornt_map_tuple=(None, )): self.weights = numpy.asarray(weights) self._intrinsic_orientation_permutation_map_tuple = io_ornt_map_tuple + def __repr__(self): + return f"{type(self).__name__}({self.point_set!r}, {self.weights!r}, {self.ref_el!r}, {self._intrinsic_orientation_permutation_map_tuple!r})" + @cached_property def point_set(self): pass # set at initialisation @@ -131,6 +147,9 @@ def __init__(self, factors, ref_el=None): for m in factor._intrinsic_orientation_permutation_map_tuple ) + def __repr__(self): + return f"{type(self).__name__}({self.factors!r}, {self.ref_el!r})" + @cached_property def point_set(self): return TensorPointSet(q.point_set for q in self.factors) diff --git a/test/finat/test_create_fiat_element.py b/test/finat/test_create_fiat_element.py index f99053ab..6530c824 100644 --- a/test/finat/test_create_fiat_element.py +++ b/test/finat/test_create_fiat_element.py @@ -59,6 +59,11 @@ def tensor_name(request): ids=lambda x: x.cellname(), scope="module") def ufl_A(request, tensor_name): + if request.param == ufl.quadrilateral: + if tensor_name == "DG": + tensor_name = "DQ" + elif tensor_name == "DG L2": + tensor_name = "DQ L2" return finat.ufl.FiniteElement(tensor_name, request.param, 1) diff --git a/test/finat/test_quadrature.py b/test/finat/test_quadrature.py new file mode 100644 index 00000000..4d1acb9f --- /dev/null +++ b/test/finat/test_quadrature.py @@ -0,0 +1,19 @@ +import pytest + +from FIAT import ufc_cell +from finat.quadrature import make_quadrature + + +@pytest.mark.parametrize( + "cell_name", + ["interval", "triangle", "interval * interval", "triangle * interval"] +) +def test_quadrature_rules_are_hashable(cell_name): + ref_cell = ufc_cell(cell_name) + quadrature1 = make_quadrature(ref_cell, 3) + quadrature2 = make_quadrature(ref_cell, 3) + + assert quadrature1 is not quadrature2 + assert hash(quadrature1) == hash(quadrature2) + assert repr(quadrature1) == repr(quadrature2) + assert quadrature1 == quadrature2 From 8e7fe86ea000eb91453caa5bec62a5719a4f97ab Mon Sep 17 00:00:00 2001 From: Connor Ward Date: Fri, 7 Feb 2025 13:50:54 +0000 Subject: [PATCH 2/5] Add safe_repr function to handle array types --- FIAT/reference_element.py | 3 ++- finat/point_set.py | 5 ++-- finat/quadrature.py | 15 ++++++++--- gem/utils.py | 54 +++++++++++++++++++++++++++++++++++++++ 4 files changed, 70 insertions(+), 7 deletions(-) diff --git a/FIAT/reference_element.py b/FIAT/reference_element.py index af5c3bf1..2cad80bf 100644 --- a/FIAT/reference_element.py +++ b/FIAT/reference_element.py @@ -26,6 +26,7 @@ from math import factorial import numpy +from gem.utils import safe_repr from recursivenodes.nodes import _decode_family, _recursive from FIAT.orientation_utils import ( @@ -185,7 +186,7 @@ def __init__(self, shape, vertices, topology): self._split_cache = {} def __repr__(self): - return f"{type(self).__name__}({self.shape!r}, {self.vertices!r}, {self.topology!r})" + return f"{type(self).__name__}({self.shape!r}, {safe_repr(self.vertices)}, {self.topology!r})" def _key(self): """Hashable object key data (excluding type).""" diff --git a/finat/point_set.py b/finat/point_set.py index 5bb7f69e..2a437c9c 100644 --- a/finat/point_set.py +++ b/finat/point_set.py @@ -1,11 +1,12 @@ import abc import hashlib +from functools import cached_property from itertools import chain, product import numpy import gem -from gem.utils import cached_property +from gem.utils import safe_repr class AbstractPointSet(abc.ABC): @@ -64,7 +65,7 @@ def __init__(self, point): self.point = point def __repr__(self): - return f"{type(self).__name__}({self.point!r})" + return f"{type(self).__name__}({safe_repr(self.point)})" @cached_property def points(self): diff --git a/finat/quadrature.py b/finat/quadrature.py index 6ea800a1..fddea657 100644 --- a/finat/quadrature.py +++ b/finat/quadrature.py @@ -1,13 +1,13 @@ import hashlib from abc import ABCMeta, abstractmethod -from functools import reduce +from functools import cached_property, reduce import gem import numpy from FIAT.quadrature import GaussLegendreQuadratureLineRule from FIAT.quadrature_schemes import create_quadrature as fiat_scheme from FIAT.reference_element import LINE, QUADRILATERAL, TENSORPRODUCT -from gem.utils import cached_property +from gem.utils import safe_repr from finat.point_set import GaussLegendrePointSet, PointSet, TensorPointSet @@ -65,7 +65,7 @@ def __hash__(self): return int.from_bytes(hashlib.md5(repr(self).encode()).digest(), byteorder="big") def __eq__(self, other): - return type(other) is type(self) and hash(other) == hash(self) + return type(other) is type(self) and repr(other) == repr(self) @abstractmethod def __repr__(self): @@ -124,7 +124,14 @@ def __init__(self, point_set, weights, ref_el=None, io_ornt_map_tuple=(None, )): self._intrinsic_orientation_permutation_map_tuple = io_ornt_map_tuple def __repr__(self): - return f"{type(self).__name__}({self.point_set!r}, {self.weights!r}, {self.ref_el!r}, {self._intrinsic_orientation_permutation_map_tuple!r})" + return ( + f"{type(self).__name__}(" + f"{self.point_set!r}, " + f"{safe_repr(self.weights)}, " + f"{self.ref_el!r}, " + f"{self._intrinsic_orientation_permutation_map_tuple!r}" + ")" + ) @cached_property def point_set(self): diff --git a/gem/utils.py b/gem/utils.py index 7e855f76..fbbbf2cf 100644 --- a/gem/utils.py +++ b/gem/utils.py @@ -1,5 +1,11 @@ import collections +import functools +import numbers from functools import cached_property # noqa: F401 +from typing import Any + +import numpy as np +from ufl.constantvalue import format_float def groupby(iterable, key=None): @@ -88,3 +94,51 @@ def __exit__(self, exc_type, exc_value, traceback): assert self.state is variable._head value, variable._head = variable._head self.state = None + + +@functools.singledispatch +def safe_repr(obj: Any) -> str: + """Return a 'safe' repr for an object, accounting for floating point error. + + Parameters + ---------- + obj : + The object to produce a repr for. + + Returns + ------- + str : + A repr for the object. + + """ + raise TypeError(f"Cannot provide a safe repr for {type(obj).__name__}") + + +@safe_repr.register(str) +def _(text: str) -> str: + return text + + +@safe_repr.register(numbers.Integral) +def _(num: numbers.Integral) -> str: + return repr(num) + + +@safe_repr.register(numbers.Real) +def _(num: numbers.Real) -> str: + return format_float(num) + + +@safe_repr.register(np.ndarray) +def _(array: np.ndarray) -> str: + return f"{type(array).__name__}([{', '.join(map(safe_repr, array))}])" + + +@safe_repr.register(list) +def _(list_: list) -> str: + return f"[{', '.join(map(safe_repr, list_))}]" + + +@safe_repr.register(tuple) +def _(tuple_: tuple) -> str: + return f"({', '.join(map(safe_repr, tuple_))})" From 1199ef11c05cab463035c02bde340c648b03d4fb Mon Sep 17 00:00:00 2001 From: Connor Ward Date: Fri, 7 Feb 2025 13:51:33 +0000 Subject: [PATCH 3/5] linting --- finat/quadrature.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/finat/quadrature.py b/finat/quadrature.py index fddea657..10e4048b 100644 --- a/finat/quadrature.py +++ b/finat/quadrature.py @@ -126,10 +126,10 @@ def __init__(self, point_set, weights, ref_el=None, io_ornt_map_tuple=(None, )): def __repr__(self): return ( f"{type(self).__name__}(" - f"{self.point_set!r}, " - f"{safe_repr(self.weights)}, " - f"{self.ref_el!r}, " - f"{self._intrinsic_orientation_permutation_map_tuple!r}" + f"{self.point_set!r}, " + f"{safe_repr(self.weights)}, " + f"{self.ref_el!r}, " + f"{self._intrinsic_orientation_permutation_map_tuple!r}" ")" ) From ea1c69c4fdd0c9413657e187ec0ed0f525d07db3 Mon Sep 17 00:00:00 2001 From: Connor Ward Date: Fri, 7 Feb 2025 14:59:21 +0000 Subject: [PATCH 4/5] better precision --- gem/utils.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/gem/utils.py b/gem/utils.py index fbbbf2cf..fae2842c 100644 --- a/gem/utils.py +++ b/gem/utils.py @@ -5,7 +5,6 @@ from typing import Any import numpy as np -from ufl.constantvalue import format_float def groupby(iterable, key=None): @@ -126,7 +125,9 @@ def _(num: numbers.Integral) -> str: @safe_repr.register(numbers.Real) def _(num: numbers.Real) -> str: - return format_float(num) + # set roundoff to half the significant figures of machine epsilon + precision = np.finfo(num).precision // 2 + return "{:.{prec}}".format(num, prec=precision) @safe_repr.register(np.ndarray) From 49627deb6e3fb689d1036e219dcd4d04e5cf596a Mon Sep 17 00:00:00 2001 From: Connor Ward Date: Fri, 7 Feb 2025 15:31:18 +0000 Subject: [PATCH 5/5] fixup --- gem/utils.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/gem/utils.py b/gem/utils.py index fae2842c..a6bf97df 100644 --- a/gem/utils.py +++ b/gem/utils.py @@ -125,8 +125,8 @@ def _(num: numbers.Integral) -> str: @safe_repr.register(numbers.Real) def _(num: numbers.Real) -> str: - # set roundoff to half the significant figures of machine epsilon - precision = np.finfo(num).precision // 2 + # set roundoff to close-to-but-not-exactly machine epsilon + precision = np.finfo(num).precision - 2 return "{:.{prec}}".format(num, prec=precision)