Simplex Domain (#152)
* add Simplex Domain * unittests for Simplex * docs update --------- Co-authored-by: Dario Coscia <93731561+dario-coscia@users.noreply.github.com>
This commit is contained in:
Ben Volokh
@@ -6,7 +6,8 @@ __all__ = [
|
||||
'Intersection',
|
||||
'Exclusion',
|
||||
'Difference',
|
||||
'OperationInterface'
|
||||
'OperationInterface',
|
||||
'SimplexDomain'
|
||||
]
|
||||
|
||||
from .location import Location
|
||||
@@ -17,3 +18,4 @@ from .intersection_domain import Intersection
|
||||
from .union_domain import Union
|
||||
from .difference_domain import Difference
|
||||
from .operation_interface import OperationInterface
|
||||
from .simplex import SimplexDomain
|
||||
|
||||
225
pina/geometry/simplex.py
Normal file
225
pina/geometry/simplex.py
Normal file
@@ -0,0 +1,225 @@
|
||||
import torch
|
||||
from .location import Location
|
||||
from pina.geometry import CartesianDomain
|
||||
from pina import LabelTensor
|
||||
from ..utils import check_consistency
|
||||
|
||||
|
||||
class SimplexDomain(Location):
|
||||
"""PINA implementation of a Simplex."""
|
||||
|
||||
def __init__(self, simplex_matrix, sample_surface=False):
|
||||
"""
|
||||
:param simplex_matrix: A matrix of LabelTensor objects representing
|
||||
a vertex of the simplex (a tensor), and the coordinates of the
|
||||
point (a list of labels).
|
||||
:type simplex_matrix: list[LabelTensor]
|
||||
:param sample_surface: A variable for choosing sample strategies. If
|
||||
`sample_surface=True` only samples on the Simplex surface
|
||||
frontier are taken. If `sample_surface=False`, no such criteria
|
||||
is followed.
|
||||
:type sample_surface: bool
|
||||
.. warning::
|
||||
Sampling for dimensions greater or equal to 10 could result
|
||||
in a shrinking of the simplex, which degrades the quality
|
||||
of the samples. For dimensions higher than 10, other algorithms
|
||||
for sampling should be used.
|
||||
:Example:
|
||||
>>> spatial_domain = SimplexDomain(
|
||||
[
|
||||
LabelTensor(torch.tensor([[0, 0]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[1, 1]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[0, 2]]), labels=["x", "y"]),
|
||||
], sample_surface = True
|
||||
)
|
||||
"""
|
||||
|
||||
# check consistency of sample_surface as bool
|
||||
check_consistency(sample_surface, bool)
|
||||
self._sample_surface = sample_surface
|
||||
|
||||
# check consistency of simplex_matrix as list or tuple
|
||||
check_consistency([simplex_matrix], (list, tuple))
|
||||
|
||||
# check everything within simplex_matrix is a LabelTensor
|
||||
check_consistency(simplex_matrix, LabelTensor)
|
||||
|
||||
# check consistency of labels
|
||||
self._coordinates = simplex_matrix[0].labels
|
||||
if not all(vertex.labels == self._coordinates for vertex in simplex_matrix):
|
||||
raise ValueError(f"Labels don't match.")
|
||||
|
||||
# creating vertices matrix
|
||||
self._vertices_matrix = torch.cat(simplex_matrix)
|
||||
self._vertices_matrix.labels = self._coordinates
|
||||
|
||||
# creating basis vectors for simplex
|
||||
self._vectors_shifted = (
|
||||
(self._vertices_matrix.T)[:, :-1] - (self._vertices_matrix.T)[:, None, -1]
|
||||
)
|
||||
|
||||
# build cartesian_bound
|
||||
self._cartesian_bound = self._build_cartesian(self.variables)
|
||||
|
||||
@property
|
||||
def coordinates(self):
|
||||
return self._coordinates
|
||||
|
||||
@property
|
||||
def variables(self):
|
||||
return self._vertices_matrix
|
||||
|
||||
def _build_cartesian(self, vertices):
|
||||
"""
|
||||
Build Cartesian border for Simplex domain to be used in sampling.
|
||||
:param vertex_matrix: matrix of vertices
|
||||
:type vertices: list[list]
|
||||
:return: Cartesian border for triangular domain
|
||||
:rtype: CartesianDomain
|
||||
"""
|
||||
|
||||
span_dict = {}
|
||||
|
||||
for i, coord in enumerate(self._coordinates):
|
||||
sorted_vertices = sorted(vertices, key=lambda vertex: vertex[i])
|
||||
# respective coord bounded by the lowest and highest values
|
||||
span_dict[coord] = [sorted_vertices[0][i], sorted_vertices[-1][i]]
|
||||
|
||||
return CartesianDomain(span_dict)
|
||||
|
||||
def is_inside(self, point, check_border=False):
|
||||
"""
|
||||
Check if a point is inside the simplex.
|
||||
Uses the algorithm described involving barycentric coordinates:
|
||||
https://en.wikipedia.org/wiki/Barycentric_coordinate_system.
|
||||
.. note::
|
||||
When ```'sample_surface'``` in the ```'__init()__'```
|
||||
is set to ```'True'```, then the method only checks
|
||||
points on the surface, and not inside the domain.
|
||||
:param point: Point to be checked.
|
||||
:type point: LabelTensor
|
||||
:param check_border: Check if the point is also on the frontier
|
||||
of the simplex, default False.
|
||||
:type check_border: bool
|
||||
:return: Returning True if the point is inside, False otherwise.
|
||||
:rtype: bool
|
||||
"""
|
||||
|
||||
if not all([label in self.coordinates for label in point.labels]):
|
||||
raise ValueError(
|
||||
"Point labels different from constructor"
|
||||
f" dictionary labels. Got {point.labels},"
|
||||
f" expected {self.variables}."
|
||||
)
|
||||
|
||||
# shift point
|
||||
point_shift = point.T - (self._vertices_matrix.T)[:, None, -1]
|
||||
|
||||
# compute barycentric coordinates
|
||||
lambda_ = torch.linalg.solve(self._vectors_shifted * 1.0, point_shift * 1.0)
|
||||
lambda_1 = 1.0 - torch.sum(lambda_)
|
||||
lambdas = torch.vstack([lambda_, lambda_1])
|
||||
|
||||
# perform checks
|
||||
if not check_border:
|
||||
return all(torch.gt(lambdas, 0.0)) and all(torch.lt(lambdas, 1.0))
|
||||
|
||||
return all(torch.ge(lambdas, 0)) and (
|
||||
any(torch.eq(lambdas, 0)) or any(torch.eq(lambdas, 1))
|
||||
)
|
||||
|
||||
def _sample_interior_randomly(self, n, variables):
|
||||
"""
|
||||
Randomly sample points inside a simplex of arbitrary
|
||||
dimension, without the boundary.
|
||||
:param int n: Number of points to sample in the shape.
|
||||
:param variables: pinn variable to be sampled, defaults to 'all'.
|
||||
:type variables: str or list[str], optional
|
||||
:return: Returns tensor of n sampled points.
|
||||
:rtype: torch.Tensor
|
||||
"""
|
||||
|
||||
# =============== For Developers ================ #
|
||||
#
|
||||
# The sampling startegy used is fairly simple.
|
||||
# First we sample a random vector from the hypercube
|
||||
# which contains the simplex. Then, if the point
|
||||
# sampled is inside the simplex, we add it as a valid
|
||||
# one.
|
||||
#
|
||||
# =============================================== #
|
||||
|
||||
sampled_points = []
|
||||
while len(sampled_points) < n:
|
||||
sampled_point = self._cartesian_bound.sample(
|
||||
n=1, mode="random", variables=variables
|
||||
)
|
||||
|
||||
if self.is_inside(sampled_point, self._sample_surface):
|
||||
sampled_points.append(sampled_point)
|
||||
return torch.cat(sampled_points, dim=0)
|
||||
|
||||
def _sample_boundary_randomly(self, n):
|
||||
"""
|
||||
Randomly sample points on the boundary of a simplex
|
||||
of arbitrary dimensions.
|
||||
:param int n: Number of points to sample in the shape.
|
||||
:return: Returns tensor of n sampled points
|
||||
:rtype: torch.Tensor
|
||||
"""
|
||||
|
||||
# =============== For Developers ================ #
|
||||
#
|
||||
# The sampling startegy used is fairly simple.
|
||||
# We first sample the lambdas in [0, 1] domain,
|
||||
# we then set to zero only one lambda, and normalize.
|
||||
# Finally, we compute the matrix product between the
|
||||
# lamdas and the vertices matrix.
|
||||
#
|
||||
# =============================================== #
|
||||
|
||||
sampled_points = []
|
||||
|
||||
while len(sampled_points) < n:
|
||||
# extract number of vertices
|
||||
number_of_vertices = self._vertices_matrix.shape[1]
|
||||
# extract idx lambda to set to zero randomly
|
||||
idx_lambda = torch.randint(low=0, high=number_of_vertices, size=(1,))
|
||||
# build lambda vector
|
||||
# 1. sampling [1, 2)
|
||||
lambdas = torch.rand((number_of_vertices, 1))
|
||||
# 2. setting lambdas[idx_lambda] to 0
|
||||
lambdas[idx_lambda] = 0
|
||||
# 3. normalize
|
||||
lambdas /= lambdas.sum()
|
||||
# 4. compute dot product
|
||||
sampled_points.append(self._vertices_matrix @ lambdas)
|
||||
|
||||
return torch.cat(sampled_points, dim=1).T
|
||||
|
||||
def sample(self, n, mode="random", variables="all"):
|
||||
"""
|
||||
Sample n points from Simplex domain.
|
||||
:param int n: Number of points to sample in the shape.
|
||||
:param str mode: Mode for sampling, defaults to 'random'.
|
||||
Available modes include: 'random'.
|
||||
:param variables: pinn variable to be sampled, defaults to 'all'.
|
||||
:type variables: str or list[str], optional
|
||||
:return: Returns LabelTensor of n sampled points
|
||||
:rtype: LabelTensor(tensor)
|
||||
.. warning::
|
||||
When ``sample_surface = True`` in the initialization, all
|
||||
the variables are sampled, despite passing different once
|
||||
in ``variables``.
|
||||
"""
|
||||
|
||||
if mode in ["random"]:
|
||||
if self._sample_surface:
|
||||
sample_pts = self._sample_boundary_randomly(n)
|
||||
else:
|
||||
sample_pts = self._sample_interior_randomly(n, variables)
|
||||
|
||||
else:
|
||||
raise NotImplementedError(f"mode={mode} is not implemented.")
|
||||
|
||||
return LabelTensor(sample_pts, labels=self.variables)
|
||||
141
tests/test_geometry/test_simplex.py
Normal file
141
tests/test_geometry/test_simplex.py
Normal file
@@ -0,0 +1,141 @@
|
||||
import torch
|
||||
import pytest
|
||||
|
||||
from pina import LabelTensor
|
||||
from pina.geometry import SimplexDomain
|
||||
|
||||
|
||||
def test_constructor():
|
||||
SimplexDomain(
|
||||
[
|
||||
LabelTensor(torch.tensor([[0, 0]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[1, 1]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[0, 2]]), labels=["x", "y"]),
|
||||
]
|
||||
)
|
||||
SimplexDomain(
|
||||
[
|
||||
LabelTensor(torch.tensor([[0, 0]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[1, 1]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[0, 2]]), labels=["x", "y"]),
|
||||
],
|
||||
sample_surface=True,
|
||||
)
|
||||
with pytest.raises(ValueError):
|
||||
SimplexDomain(
|
||||
[
|
||||
LabelTensor(torch.tensor([[0, 0]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[1, 1]]), labels=["x", "z"]),
|
||||
LabelTensor(torch.tensor([[0, 2]]), labels=["x", "a"]),
|
||||
]
|
||||
)
|
||||
SimplexDomain(
|
||||
[
|
||||
LabelTensor(torch.tensor([[0, 0]]), labels=["x", "y"]),
|
||||
[1, 1],
|
||||
LabelTensor(torch.tensor([[0, 2]]), labels=["x", "y"]),
|
||||
]
|
||||
)
|
||||
|
||||
|
||||
def test_is_inside_faulty_point():
|
||||
domain = SimplexDomain(
|
||||
[
|
||||
LabelTensor(torch.tensor([[0, 0]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[2, 2]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[4, 0]]), labels=["x", "y"]),
|
||||
],
|
||||
sample_surface=True,
|
||||
)
|
||||
pt = LabelTensor(torch.tensor([[2, 1]]), ["x", "z"])
|
||||
with pytest.raises(ValueError):
|
||||
assert domain.is_inside(point=pt, check_border=False) == True
|
||||
|
||||
|
||||
def test_is_inside_2D_check_border_true():
|
||||
domain = SimplexDomain(
|
||||
[
|
||||
LabelTensor(torch.tensor([[0, 0]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[2, 2]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[4, 0]]), labels=["x", "y"]),
|
||||
],
|
||||
sample_surface=True,
|
||||
)
|
||||
pt1 = LabelTensor(torch.tensor([[0, 0]]), ["x", "y"])
|
||||
pt2 = LabelTensor(torch.tensor([[2, 2]]), ["x", "y"])
|
||||
pt3 = LabelTensor(torch.tensor([[4, 0]]), ["x", "y"])
|
||||
pt4 = LabelTensor(torch.tensor([[3, 1]]), ["x", "y"])
|
||||
pt5 = LabelTensor(torch.tensor([[2, 1]]), ["x", "y"])
|
||||
pt6 = LabelTensor(torch.tensor([[100, 100]]), ["x", "y"])
|
||||
pts = [pt1, pt2, pt3, pt4, pt5, pt6]
|
||||
for pt, exp_result in zip(pts, [True, True, True, True, False, False]):
|
||||
assert domain.is_inside(point=pt, check_border=True) == exp_result
|
||||
|
||||
|
||||
def test_is_inside_2D_check_border_false():
|
||||
domain = SimplexDomain(
|
||||
[
|
||||
LabelTensor(torch.tensor([[0, 0]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[2, 2]]), labels=["x", "y"]),
|
||||
LabelTensor(torch.tensor([[4, 0]]), labels=["x", "y"]),
|
||||
],
|
||||
sample_surface=False,
|
||||
)
|
||||
pt1 = LabelTensor(torch.tensor([[0, 0]]), ["x", "y"])
|
||||
pt2 = LabelTensor(torch.tensor([[2, 2]]), ["x", "y"])
|
||||
pt3 = LabelTensor(torch.tensor([[4, 0]]), ["x", "y"])
|
||||
pt4 = LabelTensor(torch.tensor([[3, 1]]), ["x", "y"])
|
||||
pt5 = LabelTensor(torch.tensor([[2, 1]]), ["x", "y"])
|
||||
pt6 = LabelTensor(torch.tensor([[2.5, 1]]), ["x", "y"])
|
||||
pt7 = LabelTensor(torch.tensor([[100, 100]]), ["x", "y"])
|
||||
pts = [pt1, pt2, pt3, pt4, pt5, pt6, pt7]
|
||||
for pt, exp_result in zip(pts, [False, False, False, False, True, True, False]):
|
||||
assert domain.is_inside(point=pt, check_border=False) == exp_result
|
||||
|
||||
|
||||
def test_is_inside_3D_check_border_true():
|
||||
domain = SimplexDomain(
|
||||
[
|
||||
LabelTensor(torch.tensor([[0, 0, 0]]), labels=["x", "y", "z"]),
|
||||
LabelTensor(torch.tensor([[2, 2, 0]]), labels=["x", "y", "z"]),
|
||||
LabelTensor(torch.tensor([[4, 0, 0]]), labels=["x", "y", "z"]),
|
||||
LabelTensor(torch.tensor([[0, 0, 20]]), labels=["x", "y", "z"]),
|
||||
],
|
||||
sample_surface=True,
|
||||
)
|
||||
pt1 = LabelTensor(torch.tensor([[0, 0, 0]]), ["x", "y", "z"])
|
||||
pt2 = LabelTensor(torch.tensor([[2, 2, 0]]), ["x", "y", "z"])
|
||||
pt3 = LabelTensor(torch.tensor([[4, 0, 0]]), ["x", "y", "z"])
|
||||
pt4 = LabelTensor(torch.tensor([[3, 1, 0]]), ["x", "y", "z"])
|
||||
pt5 = LabelTensor(torch.tensor([[2, 1, 0]]), ["x", "y", "z"])
|
||||
pt6 = LabelTensor(torch.tensor([[100, 100, 0]]), ["x", "y", "z"])
|
||||
pt7 = LabelTensor(torch.tensor([[0, 0, 19]]), ["x", "y", "z"])
|
||||
pt8 = LabelTensor(torch.tensor([[0, 0, 20]]), ["x", "y", "z"])
|
||||
pt9 = LabelTensor(torch.tensor([[2, 1, 1]]), ["x", "y", "z"])
|
||||
pts = [pt1, pt2, pt3, pt4, pt5, pt6, pt7, pt8, pt9]
|
||||
for pt, exp_result in zip(
|
||||
pts, [True, True, True, True, True, False, True, True, False]
|
||||
):
|
||||
assert domain.is_inside(point=pt, check_border=True) == exp_result
|
||||
|
||||
|
||||
def test_is_inside_3D_check_border_false():
|
||||
domain = SimplexDomain(
|
||||
[
|
||||
LabelTensor(torch.tensor([[0, 0, 0]]), labels=["x", "y", "z"]),
|
||||
LabelTensor(torch.tensor([[2, 2, 0]]), labels=["x", "y", "z"]),
|
||||
LabelTensor(torch.tensor([[4, 0, 0]]), labels=["x", "y", "z"]),
|
||||
LabelTensor(torch.tensor([[0, 0, 20]]), labels=["x", "y", "z"]),
|
||||
],
|
||||
sample_surface=False,
|
||||
)
|
||||
pt1 = LabelTensor(torch.tensor([[0, 0, 0]]), ["x", "y", "z"])
|
||||
pt2 = LabelTensor(torch.tensor([[3, 1, 0]]), ["x", "y", "z"])
|
||||
pt3 = LabelTensor(torch.tensor([[2, 1, 0]]), ["x", "y", "z"])
|
||||
pt4 = LabelTensor(torch.tensor([[100, 100, 0]]), ["x", "y", "z"])
|
||||
pt5 = LabelTensor(torch.tensor([[0, 0, 19]]), ["x", "y", "z"])
|
||||
pt6 = LabelTensor(torch.tensor([[0, 0, 20]]), ["x", "y", "z"])
|
||||
pt7 = LabelTensor(torch.tensor([[2, 1, 1]]), ["x", "y", "z"])
|
||||
pts = [pt1, pt2, pt3, pt4, pt5, pt6, pt7]
|
||||
for pt, exp_result in zip(pts, [False, False, False, False, False, False, True]):
|
||||
assert domain.is_inside(point=pt, check_border=False) == exp_result
|
||||
Reference in New Issue
Block a user