from collections.abc import Container, Sequence # , Iterable
import numpy as np
from ddg.nets._traverser import Traverser
def _is_bounded_interval(interval):
"""Checks whether an interval is bounded."""
return interval[0] > -np.inf and interval[1] < np.inf
def _domain_decomposer(domain):
"""
Decomposes a list of lists describing a domain into its intervals
and periodicity.
"""
# TODO: handle cases of a parametrizable domains given for example as
# [['a',0,True],['b','c']]
intervals = [x[:2] for x in domain]
periodicity = [x[-1] if isinstance(x[-1], bool) else False for x in domain]
return intervals, periodicity
[docs]class Domain(Container):
"""
Base class for storing information on a `point-domain` of a net.
Parameters
----------
dimension : non negative int
Dimension of the Domain.
Attributes
----------
dimension : non negative int
Dimension of the Domain.
"""
def __init__(self, dimension, *args):
self.dimension = dimension
def __contains__(self, point):
raise NotImplementedError
[docs]class EmptyDomain(Domain):
def __init__(self):
super().__init__(None)
def __contains__(self, point):
return False
[docs]class SmoothDomain(Domain):
"""
Factory class for the smooth domains.
Parameters
----------
domain_type : {'rectangular'}, optional
The type of discrete domain to use.
Raises
------
NotImplementedError
If the domain type does not exist.
See Also
--------
Domain
SmoothRectangularDomain
"""
def __new__(cls, *args, domain_type="rectangular", **kwargs):
"""Create a smooth domain of the given type."""
# TODO: Better handling of optional arguments!!!
if domain_type == "rectangular":
return SmoothRectangularDomain.__new__(cls, *args, **kwargs)
else:
raise NotImplementedError(
"There is no domain type named: {}".format(domain_type)
)
[docs]class SmoothRectangularDomain(SmoothDomain, Sequence):
"""
Store information on a rectangular domain.
Parameters
----------
domain : list of tuples (float, float[, bool])
Tuples define the intervals as (lowerbound, upperbound, periodic)
Periodic is optional and assumed False if not given.
Attributes
----------
periodicity : set
Set containing indices of periodic coordinate directions.
intervals : list of tuples
List of the coordinate intervals.
bounded_directions : list
List containing the indices of the bounded directions of the domain.
unbounded_directions : list
List containing the indices of the unbounded directions of the domain.
Methods
-------
recover
Returns the current state of the domain in a format that can be used
to initialize a constant domain with the same intervals and periodicity
See Also
--------
SmoothDomain
Domain
"""
def __new__(cls, domain):
if len(domain) == 1:
return super(SmoothDomain, cls).__new__(SmoothInterval)
else:
return super(SmoothDomain, cls).__new__(SmoothRectangularDomain)
def __init__(self, domain):
intervals, periodicity = _domain_decomposer(domain)
self._intervals = [[float(el) for el in interval] for interval in intervals]
self._periodicity = set(i for i, val in enumerate(periodicity) if val)
def __getitem__(self, key):
return self.intervals[key]
def __len__(self):
return len(self.intervals)
def __contains__(self, point):
for [a, b], x in zip(self.intervals, point):
if not (a <= x <= b):
return False
return True
def __repr__(self):
return self.__class__.__name__ + "(" + repr(self.recover()) + ")"
@property
def intervals(self):
return self._intervals
@property
def periodicity(self):
return self._periodicity
@property
def bounded(self):
return all([_is_bounded_interval(x) for x in self.intervals])
@property
def dimension(self):
return len(self.intervals)
[docs] def recover(self):
dlist = [x + [i in self.periodicity] for i, x in enumerate(self.intervals)]
return dlist
@property
def unbounded_directions(self):
unbounded = []
for idx, el in enumerate(self.intervals):
if not _is_bounded_interval(el):
unbounded.append(idx)
return unbounded
@property
def bounded_directions(self):
bounded = []
for idx, el in enumerate(self.intervals):
if _is_bounded_interval(el):
bounded.append(idx)
return bounded
[docs]class SmoothInterval(SmoothRectangularDomain):
"""
Special case of a one-dimensional rectangular domain.
Parameters
----------
interval : tuple of two floats and an optional bool
A tuple defining the endpoints and periodicity of the interval.
Periodicity is assumed False if not given.
Attributes
----------
interval: tuple
The tuple defining the interval
See Also
--------
SmoothRectangularDomain
"""
def __new__(cls, interval):
return super(SmoothDomain, cls).__new__(SmoothInterval)
def __init__(self, interval):
if len(interval) == 1:
super().__init__(interval)
else:
super().__init__([interval])
@property
def interval(self):
"""Interval as a tuple."""
return self.intervals[0]
@property
def periodic(self):
return 0 in self.periodicity
[docs] def recover(self):
return self.interval + [self.periodic]
[docs]class DiscreteDomain(Domain): # , Iterable):
"""
Factory class for the discrete domains.
Parameters
----------
traverser : ddg.nets.traverser.Traverser
Traversing scheme of the discrete net.
domain_type : {'rectangular', 'triangular', 'diagonal'}, optional
The type of discrete domain to use.
Attributes
----------
edge_data : list
List of edges in the domain given by index pairs.
face_data : list
List of edges in the domain given by tuples of indices
(at least triples).
traverser : ddg.nets.traverser.Traverser
Traversing scheme of the discrete net.
Raises
------
NotImplementedError
If the domain type does not exist.
See Also
--------
Domain, DiscreteRectangularDomain, DiscreteTriangularDomain,
DiscreteDiagonalDomain
"""
def __new__(cls, *args, domain_type="rectangular", **kwargs):
"""Create a discrete domain of the given type."""
if domain_type == "rectangular":
return DiscreteRectangularDomain.__new__(cls, *args, **kwargs)
elif domain_type == "triangular":
return super().__new__(DiscreteTriangularDomain)
elif domain_type == "diagonal":
return super().__new__(DiscreteDiagonalDomain)
else:
raise NotImplementedError(
"There is no domain type named: {}".format(domain_type)
)
def __init__(self, traverser, *args):
self.traverser = traverser
self._edge_data = None
self._face_data = None
super().__init__(*args)
def __contains__(self, val):
raise NotImplementedError
@property
def edge_data(self):
if not self._edge_data:
self._gen_edge_data()
return self._edge_data
@property
def face_data(self):
if not self._face_data:
self._gen_face_data()
return self._face_data
def _gen_face_data(self):
raise NotImplementedError
def _gen_edge_data(self):
raise NotImplementedError
[docs]class DiscreteRectangularDomain(DiscreteDomain, Sequence):
"""
Store data on a rectangular index domain.
Parameters
----------
*args : tuples of int
List of ranges of vertices.
periodicity : set or list, optional
Specify which coordinate directions are periodic.
(default: No periodicity).
traverser : ddg.nets.traverser.Traverser, optional
Traversing scheme of the discrete domain.
Attributes
----------
periodicity : set
Set containing indices of periodic coordinate directions.
intervals : list of tuples
List of the coordinate intervals.
bounded_directions : list
List containing the indices of the bounded directions of the domain.
unbounded_directions : list
List containing the indices of the unbounded directions of the domain.
Methods
-------
recover
Returns the current state of the domain in a format that can be used
to initialize a constant domain with the same intervals and periodicity
See Also
--------
Domain, DiscreteDomain
"""
def __new__(cls, domain, **kwargs):
if len(domain) == 1:
return super(DiscreteDomain, cls).__new__(DiscreteInterval)
else:
return super(DiscreteDomain, cls).__new__(DiscreteRectangularDomain)
def __init__(self, domain, traverser=None):
intervals, periodicity = _domain_decomposer(domain)
self._intervals = [[el for el in interval] for interval in intervals]
self._periodicity = set(i for i, val in enumerate(periodicity) if val)
self._double_edged = False
for i in self._periodicity:
interval = self.intervals[i]
if interval[1] - interval[0] == 1:
self._double_edged = True
# temp, periodicity = domain_decomposer(domain)
# self.periodicity = set(i for i,val in enumerate(periodicity) if val)
# self.intervals = domain
# if traverser is None:
# traverser = Traverser.Linear(intervals=self.intervals)
# super().__init__(traverser, len(self.intervals))
def __getitem__(self, key):
return self.intervals[key]
def __len__(self):
return len(self.intervals)
def __contains__(self, point):
for [a, b], x in zip(self.intervals, point):
if not (a <= x <= b):
return False
return True
def __repr__(self):
return (
self.__class__.__name__
+ "("
+ ", ".join([repr(x) for x in self.intervals])
+ ")"
)
def _gen_face_data(self):
if sum([b - a + 1 for a, b in self.intervals]) == float("inf"):
raise NotImplementedError(
"Face data generation not implemented for unbounded domains."
)
if self.dimension == 1:
self._face_data = []
elif self.dimension == 2:
width, height = [b - a + 1 for a, b in self.intervals]
idx = self.traverser.idx
shift = (self.intervals[0][0], self.intervals[1][0])
self._face_data = []
for k in range(width - 1):
for l in range(height - 1):
i, j = k + shift[0], l + shift[1]
self._face_data.append(
(idx(i, j), idx(i + 1, j), idx(i + 1, j + 1), idx(i, j + 1))
)
for p in self.periodicity:
p_start, p_end = self.intervals[p]
start, end = self.intervals[(p + 1) % 2]
for i in range(start, end):
tmp_idx1, tmp_idx2 = [0, 0], [0, 0]
tmp_idx3, tmp_idx4 = [0, 0], [0, 0]
tmp_idx1[p], tmp_idx1[(p + 1) % 2] = p_start, i
tmp_idx2[p], tmp_idx2[(p + 1) % 2] = p_end, i
tmp_idx3[p], tmp_idx3[(p + 1) % 2] = p_end, i + 1
tmp_idx4[p], tmp_idx4[(p + 1) % 2] = p_start, i + 1
if p == 1:
self._face_data.append(
(
idx(*tmp_idx1),
idx(*tmp_idx2),
idx(*tmp_idx3),
idx(*tmp_idx4),
)
)
else:
self._face_data.append(
(
idx(*tmp_idx1),
idx(*tmp_idx4),
idx(*tmp_idx3),
idx(*tmp_idx2),
)
)
if len(self.periodicity) == 2:
(a, b), (c, d) = self.intervals
self._face_data.append((idx(b, d), idx(a, d), idx(a, c), idx(b, c)))
else:
idx = self.traverser.idx
lower = [self.intervals[i][0] for i in range(self.dimension)]
upper = np.array([self.intervals[i][-1] for i in range(self.dimension)])
self._face_data = []
for point in self.traverser:
pidx = idx(*point)
mask = upper == point
for i in range(len(mask)):
for j in range(i + 1, len(mask)):
t1 = None
t2 = None
t3 = None
if not (mask[i] or mask[j]):
t1 = np.array(point)
t2 = np.array(point)
t3 = np.array(point)
t1[i] += 1
t3[i] += 1
t2[j] += 1
t3[j] += 1
elif mask[i] and mask[j]:
if i in self.periodicity and j in self.periodicity:
t1 = np.array(point)
t2 = np.array(point)
t3 = np.array(point)
t1[i] = lower[i]
t3[i] = lower[i]
t2[j] = lower[j]
t3[j] = lower[j]
elif mask[i] and i in self.periodicity:
t1 = np.array(point)
t2 = np.array(point)
t3 = np.array(point)
t1[i] = lower[i]
t3[i] = lower[i]
t2[j] += 1
t3[j] += 1
elif mask[j] and j in self.periodicity:
t1 = np.array(point)
t2 = np.array(point)
t3 = np.array(point)
t1[i] += 1
t3[i] += 1
t2[j] = lower[j]
t3[j] = lower[j]
if t1 is not None:
self._face_data.append((pidx, idx(*t1), idx(*t3), idx(*t2)))
def _gen_edge_data(self):
if sum([b - a + 1 for a, b in self.intervals]) == float("inf"):
raise NotImplementedError
if self.dimension == 1:
length = self.intervals[0][1] - self.intervals[0][0] + 1
idx = self.traverser.idx
shift = self.intervals[0][0]
self._edge_data = []
if 0 in self.periodicity:
end = shift + length
else:
end = shift + length - 1
for i in range(shift, end):
self._edge_data.append((idx(i), idx(i + 1)))
elif self.dimension == 2:
width, height = [b - a + 1 for a, b in self.intervals]
idx = self.traverser.idx
shift = (self.intervals[0][0], self.intervals[1][0])
self._edge_data = []
for i in range(shift[0], shift[0] + width - 1):
for j in range(shift[1], shift[1] + height - 1):
self._edge_data.append((idx(i, j), idx(i + 1, j)))
self._edge_data.append((idx(i, j), idx(i, j + 1)))
j = shift[1] + height - 1
for i in range(shift[0], shift[0] + width - 1):
self._edge_data.append((idx(i, j), idx(i + 1, j)))
i = shift[0] + width - 1
for j in range(shift[1], shift[1] + height - 1):
self._edge_data.append((idx(i, j), idx(i, j + 1)))
for p in self.periodicity:
p_start, p_end = self.intervals[p]
start, end = self.intervals[(p + 1) % 2]
for i in range(start, end + 1):
tmp_idx1, tmp_idx2 = [0, 0], [0, 0]
tmp_idx1[p], tmp_idx1[(p + 1) % 2] = p_start, i
tmp_idx2[p], tmp_idx2[(p + 1) % 2] = p_end, i
self._edge_data.append((idx(*tmp_idx1), idx(*tmp_idx2)))
else:
idx = self.traverser.idx
lower = [self.intervals[i][0] for i in range(self.dimension)]
upper = np.array([self.intervals[i][-1] for i in range(self.dimension)])
self._edge_data = []
for point in self.traverser:
pidx = idx(*point)
mask = upper == point
for i in range(len(mask)):
temp = None
if not mask[i]:
temp = np.array(point)
temp[i] = temp[i] + 1
elif i in self.periodicity:
temp = np.array(point)
temp[i] = lower[i]
if temp is not None:
self._edge_data.append((pidx, idx(*temp)))
@property
def edge_data(self):
if not hasattr(self, "_edge_data"):
self._gen_edge_data()
return self._edge_data
@property
def face_data(self):
if not hasattr(self, "_face_data"):
self._gen_face_data()
return self._face_data
@property
def dimension(self):
return len(self.intervals)
@property
def double_edged(self):
return self._double_edged
@property
def intervals(self):
return self._intervals
@property
def periodicity(self):
return self._periodicity
# TODO: Non-linear Traversers and handling if domain has not changed.
@property
def traverser(self):
return Traverser.Linear(intervals=self.intervals)
@property
def bounded(self):
return all([_is_bounded_interval(x) for x in self.intervals])
[docs] def recover(self):
dlist = [x + [i in self.periodicity] for i, x in enumerate(self.intervals)]
return dlist
@property
def unbounded_directions(self):
unbounded = []
for idx, el in enumerate(self.intervals):
if not _is_bounded_interval(el):
unbounded.append(idx)
return unbounded
@property
def bounded_directions(self):
bounded = []
for idx, el in enumerate(self.intervals):
if _is_bounded_interval(el):
bounded.append(idx)
return bounded
[docs]class DiscreteInterval(DiscreteRectangularDomain):
"""
Special case of a one-dimensional rectangular domain.
Parameters
----------
interval : tuple of two floats [and a single bool] or
function returning such a tuple
One tuple defining the endpoints of the interval and if
it is periodic or not.
Attributes
----------
interval: tuple
The tuple defining the interval
periodic : bool
See Also
--------
DiscreteRectangularDomain
"""
def __new__(cls, *args):
return super(DiscreteDomain, cls).__new__(DiscreteInterval)
def __init__(self, interval):
if len(interval) == 1:
super().__init__(interval)
else:
super().__init__([interval])
@property
def interval(self):
"""Interval as a tuple."""
return self.intervals[0]
@property
def periodic(self):
return 0 in self.periodicity
[docs] def recover(self):
return self.interval + [self.periodic]
[docs]class DiscreteTriangularDomain(DiscreteDomain):
r"""
Store data on a triangular index domain.
Parameters
----------
base : float
Number of vertex points in base edge of triangle.
shift : tuple of int, optional
Base shift for triangle.
mode : {1, ..., 8}, optional
Select triangle shape.
Attributes
----------
edge_data : list
List of edges in the domain given by index pairs.
face_data : list
List of edges in the domain given by tuples of indices
(at least triples).
Notes
-----
* The number of vertex points is calculated from `base` by ceiling
the absolut value of `base`.
* The triangle is constructed from `base` by::
o ^
/ | .
/ | . base+1
/ | .
o-------o v
< . . . >
base+1
* If `base` is `inf` the half axis is taken as the base edge.
* The `mode` selects one of the following triangles (NOT IMPLEMENTED YET)::
o--------o--------o
| \ 3 | 2 / |
| \ | / |
| 4 \ | / 1 |
o--------o--------o
| 5 / | \ 8 |
| / | \ |
| / 6 | 7 \ |
o--------o--------o
"""
def __init__(self, base, *args, shift=(0, 0), mode=1, **kwargs):
# TODO: Use `mode` option.
self._base = np.ceil(np.abs(base))
self._shift = shift
super().__init__(
Traverser.Linear(self._base, shape="triang", dirct="S", shift=self._shift),
2,
)
self._gen_face_data()
self._gen_edge_data()
def __contains__(self, idx):
width_const = self._shift[0] <= idx[0] <= self._shift[0] + self._base
height_const = self._shift[1] <= idx[1] <= self._shift[1] + idx[0]
return width_const and height_const
def _gen_face_data(self):
if self._base == float("inf"):
raise NotImplementedError
else:
idx = self.traverser.idx
self._base = int(self._base)
self._face_data = []
for k in range(1, self._base - 1):
for l in range(k):
i, j = k + self._shift[0], l + self._shift[1]
self._face_data.append(
(idx(i, j), idx(i + 1, j), idx(i + 1, j + 1), idx(i, j + 1))
)
# add triangles (not quads) at diagonal
# TODO: Check if this is a good idea.
for k in range(self._base - 1):
i, j = k + self._shift[0], k + self._shift[1]
self._face_data.append((idx(i, j), idx(i + 1, j), idx(i + 1, j + 1)))
def _gen_edge_data(self):
if self._base == float("inf"):
raise NotImplementedError
else:
idx = self.traverser.idx
self._base = int(self._base)
self._edge_data = []
for k in range(1, self._base - 1):
for l in range(k):
i, j = k + self._shift[0], l + self._shift[1]
self._edge_data.append((idx(i, j), idx(i + 1, j)))
self._edge_data.append((idx(i, j), idx(i, j + 1)))
for k in range(self._base - 1):
i, j = k + self._shift[0], k + self._shift[1]
self._edge_data.append((idx(i, j), idx(i + 1, j)))
# add diagonal
# TODO: Check if this is a good idea.
self._edge_data.append((idx(i, j), idx(i + 1, j + 1)))
for k in range(self._base - 1):
i, j = self._base - 1 + self._shift[0], k + self._shift[1]
self._edge_data.append((idx(i, j), idx(i, j + 1)))
[docs]class DiscreteDiagonalDomain(DiscreteDomain):
r"""
Implement a diagonal strip index domain.
Parameters
----------
height : float
Hight of the strip.
width : float
Width of the strip.
shift : tuple of int, optional
Base shift for diagonal domain.
Attributes
----------
edge_data : list
List of edges in the domain given by index pairs.
face_data : list
List of edges in the domain given by tuples of indices
(at least triples).
Notes
-----
* The number of vertex points is calculated from the elements of
`args` by ceiling their absolut values, respectively.
* There needs to be at least one positional argument and at most
two. If `args` has length one both sides are given that lenth.
* The elements of `args` are assigned to the sides of the diagonal
strip as follows::
^ o------o
args[0]+1 . / /
. / /
v o------o
< . . .>
args[1]+1
* If `base` is `inf` the half axis is taken as the base edge.
* The `mode` selects one of the following diagonal strips
(NOT IMPLEMENTED YET)::
o------o------o o------o------o
\ \ \ / / /
\ 4 \ 3 \ / 2 / 1 /
\ \ \ / / /
o------o------o o------o------o
"""
def __init__(self, height, width, *args, mode=1, shift=(0, 0), **kwargs):
self._height, self._width = height, width
self._shift = shift
if self._width == float("inf") and self._height == float("inf"):
super().__init__(Traverser.Circular(shift=self._shift), 2)
elif self._width == float("inf") or self._height == float("inf"):
# TODO: Delete this case when `Traverser.Diag` implemented
# compatibility with infinite slide lengths.
raise NotImplementedError
else:
super().__init__(
Traverser.Diag(
self._height, self._width, shape="half", shift=self._shift
),
2,
)
def __contains__(self, idx):
i, j = idx[0] - self._shift[0], idx[1] - self._shift[1]
width_const = j <= i <= self._width + j
height_const = 0 <= j <= self._height
return width_const and height_const
def _gen_face_data(self):
if self._width == float("inf") or self._height == float("inf"):
raise NotImplementedError
else:
idx = self.traverser.idx
self._width, self._height = int(self._width), int(self._height)
self._face_data = []
for k in range(1, self._width - 1):
for l in range(0, self._height - 1):
i, j = k + l + self._shift[0], l + self._shift[1]
self._face_data.append(
(idx(i, j), idx(i + 1, j), idx(i + 1, j + 1), idx(i, j + 1))
)
# add triangles (not quads) at diagonals
# TODO: Check if this is a good idea.
for k in range(self._height - 1):
# left diagonal
i, j = k + self._shift[0], k + self._shift[1]
self._face_data.append((idx(i, j), idx(i + 1, j), idx(i + 1, j + 1)))
# right diagonal
i, j = k + self._width + self._shift[0] - 1, k + self._shift[1]
self._face_data.append((idx(i, j), idx(i + 1, j + 1), idx(i, j + 1)))
def _gen_edge_data(self):
if self._width == float("inf") or self._height == float("inf"):
raise NotImplementedError
else:
idx = self.traverser.idx
self._width, self._height = int(self._width), int(self._height)
self._edge_data = []
for k in range(1, self._width - 1):
for l in range(0, self._height - 1):
i, j = k + l + self._shift[0], l + self._shift[1]
self._edge_data.append((idx(i, j), idx(i + 1, j)))
self._edge_data.append((idx(i, j), idx(i, j + 1)))
for k in range(self._height - 1):
i, j = k + self._shift[0], k + self._shift[1]
self._edge_data.append((idx(i, j), idx(i + 1, j)))
# add left diagonal
# TODO: Check if this is a good idea.
self._edge_data.append((idx(i, j), idx(i + 1, j + 1)))
i, j = k + self._width - 1 + self._shift[0], k + self._shift[1]
self._edge_data.append((idx(i, j), idx(i, j + 1)))
# add right diagonal
# TODO: Check if this is a good idea.
self._edge_data.append((idx(i, j), idx(i + 1, j + 1)))
for k in range(self._width - 1):
i, j = (
k + self._height - 1 + self._shift[0],
self._height - 1 + self._shift[1],
)
self._edge_data.append((idx(i, j), idx(i + 1, j)))