import itertools
from collections.abc import Iterable
from functools import partial, wraps
import numpy as np
import ddg
import ddg.datastructures.nets.conversion as netsconv
import ddg.math.projective as putils
from ddg import nonexact
from ddg.datastructures.nets.domain import (
DiscreteDomain,
DiscreteRectangularDomain,
SmoothDomain,
SmoothRectangularDomain,
is_bounded_interval,
)
from ddg.datastructures.nets.net import (
DiscreteCurve,
DiscreteNet,
EmptyNet,
Net,
NetCollection,
PointNet,
SmoothCurve,
SmoothNet,
)
from ddg.math import euclidean
[docs]def domain_util(f):
"""
Decorator for domain utilities.
It is assumed that the decorated function returns a list of intervals and
that its first argument is the domain to work on.
The newly constructed domain will be of the same type as the input and for
nets and net collections the wrapped function will change their domain
inplace, i.e. replace the domain of the net with the newly constructed one.
Parameters
----------
f : Callable
domain utility function
Returns
-------
Callable : wrapped function
See Also
--------
net_util
"""
@wraps(f)
def wrapped(domain, *args, **kwargs):
if not isinstance(domain, NetCollection):
net = None
if isinstance(domain, ddg.nets.net.Net):
net = domain
domain = net.domain
# Remove cases once Smooth/DiscreteInterval is removed
if isinstance(domain, SmoothDomain):
target_type = SmoothDomain
elif isinstance(domain, DiscreteDomain):
target_type = DiscreteDomain
else:
target_type = domain.__class__
if isinstance(domain, (DiscreteDomain, SmoothDomain)):
new_domain = target_type(f(domain, *args, **kwargs))
else:
new_domain = target_type()
if net is None:
return new_domain
else:
net.domain = new_domain
return net
else:
for net in domain._nets:
wrapped(net, *args, **kwargs)
return domain
return wrapped
[docs]def net_util(f):
"""
Decorator for net utilities.
It is assumed that the decorated function returns a tuple of a function and
a list of intervals.
The newly constructed net will be if the same type as the input.
Parameters
----------
f : Callable
net utility function
Returns
-------
Callable : wrapped function
See Also
--------
domain_util
"""
@wraps(f)
def wrapped(net, *args, **kwargs):
if not isinstance(net, NetCollection):
# Remove cases once Smooth/DiscreteCurve is removed
if isinstance(net, SmoothNet):
target_type_domain = SmoothDomain
target_type_net = SmoothNet
elif isinstance(net, DiscreteNet):
target_type_domain = DiscreteDomain
target_type_net = DiscreteNet
else:
target_type_net = net.__class__
if not isinstance(net, (EmptyNet, PointNet)):
new_domain = target_type_domain(f(net, *args, **kwargs)[1])
new_net = target_type_net(f(net, *args, **kwargs)[0], new_domain)
else:
new_net = target_type_net(f(net, *args, **kwargs)[0])
return new_net
else:
new_nets = []
for _net in net._nets:
# Remove cases once Smooth/DiscreteCurve is removed
if isinstance(_net, SmoothNet):
target_type_domain = SmoothDomain
target_type_net = SmoothNet
elif isinstance(_net, DiscreteNet):
target_type_domain = DiscreteDomain
target_type_net = DiscreteNet
else:
target_type_net = _net.__class__
if not isinstance(net, (PointNet, EmptyNet)):
new_domain = target_type_domain(f(_net, *args, **kwargs)[1])
new_nets.append(
target_type_net(f(_net, *args, **kwargs)[0], new_domain)
)
else:
new_nets.append(target_type_net(f(_net, *args, **kwargs)[0]))
return NetCollection(new_nets)
return wrapped
[docs]def compose(f, n):
"""The composition x -> f(n(x)) adapted to nets and net collections.
Parameters
----------
f : callable or numpy.ndarray
n : Net or NetCollection
Returns
-------
Net or NetCollection
The same type as `n`.
Examples
--------
>>> import numpy as np
>>> from ddg.datastructures.nets.net import SmoothNet
>>> from ddg.datastructures.nets.utils import compose
>>> line = SmoothNet(lambda t: np.array([t, 0]), [(0, np.pi, True)])
>>> line(2)
array([2, 0])
>>> translated_line = compose(lambda x: x + (0, 1), line)
>>> translated_line(2)
array([2, 1])
>>> linearly_transformed_line = compose(np.array([[2, 1], [0, 1]]), line)
>>> linearly_transformed_line(2)
array([4, 0])
Unlike
.. code-block:: python
def f_after_n(x):
return f(n(x))
``compose`` retains the net domain(s) by returning a net or net collection.
>>> type(linearly_transformed_line)
<class 'ddg.datastructures.nets.net.SmoothCurve'>
>>> linearly_transformed_line.domain
SmoothInterval([0.0, 3.14..., True])
The parameter `n` is required to be a net or a net collection. If it isn't,
then there is no need to retain a domain, so
.. code-block:: python
def f_after_n(x):
return f(n(x))
suffices.
See Also
--------
ddg.datastructures.nets.utils.embed
ddg.datastructures.nets.utils.homogenize
ddg.datastructures.nets.utils.dehomogenize
"""
if isinstance(f, np.ndarray):
def g(*args, **kwargs):
return f @ n.fct(*args, **kwargs)
elif callable(f):
def g(*args, **kwargs):
return f(n.fct(*args, **kwargs))
else:
raise TypeError("f must be callable or an instance of numpy.ndarray")
# NetCollection isn't actually a Net.
# EmptyNet and PointNet have constructors that require
# points/coordinates/values rather than functions.
if isinstance(n, NetCollection):
return NetCollection([compose(f, subnet) for subnet in n], n.name)
elif isinstance(n, EmptyNet):
return EmptyNet(g(), n.name)
elif isinstance(n, PointNet):
return PointNet(g(), n.name)
elif isinstance(n, Net):
return type(n)(g, n.domain, name=n.name)
else:
raise TypeError(f"n must be an instance of {Net}")
[docs]def homogenize(net):
"""
Homogenize the values of given net.
Parameters
----------
net : ddg.datastructures.nets.net.Net
The net to operate on.
Returns
-------
A net of the same type.
See Also
--------
Net
Notes
-----
* It is assumed that the values of the net are of type
`np.array`.
"""
return compose(putils.homogenize, net)
[docs]def dehomogenize(net):
"""
Dehomogenize the values of given net.
Parameters
----------
net : ddg.datastructures.nets.net.Net
The net to operate on.
Returns
-------
A net of the same type.
See Also
--------
Net
Notes
-----
* It is assumed that the values of the net are of type
`np.array`.
"""
return compose(putils.dehomogenize, net)
[docs]def embed(net, level=0, component=-1):
"""
Embeds the given net in a space one dimension higher than before.
Parameters
----------
net : ddg.datastructures.nets.net.Net
The net to operate on.
level : float
(default 0)
component : int (default=-1)
Component to add
Returns
-------
A net of the same type.
See Also
--------
Net
Notes
-----
* It is assumed that the values of the net are of type
`np.array`.
"""
return compose(partial(euclidean.embed, index=component, level=level), net)
[docs]def vertices(net):
"""
Get all vertex values from a given net or net collection.
Parameters
----------
net : ddg.datastructures.nets.net.DiscreteNet or NetCollection
The net or net collection to obtain all the vertices from
Returns
-------
Generator
Generator containing all vertex values.
"""
# Check whether domain of net is bounded
if isinstance(net, DiscreteNet) and not net.domain.bounded:
raise TypeError("Domain of net is unbounded.")
elif isinstance(net, NetCollection):
for n in net:
if isinstance(n, DiscreteNet) and not n.domain.bounded:
raise TypeError("Domain of net is unbounded.")
if (
isinstance(net, DiscreteNet)
or isinstance(net, PointNet)
or isinstance(net, EmptyNet)
):
return iter(net)
elif isinstance(net, NetCollection):
return itertools.chain(*map(vertices, net))
else:
raise TypeError("Only supports DiscretNet and NetCollection.")
[docs]def applicable_to_netcollection(f, creates_new_nets=True):
if creates_new_nets is True:
@wraps(f)
def wrapped_f(nets, *args, **kwargs):
if isinstance(nets, NetCollection):
if kwargs.get("name") is not None:
nc = NetCollection(name=kwargs["name"])
else:
nc = NetCollection()
for net in nets:
nc.add(f(net, *args, **kwargs))
if nc._nets == nets._nets:
return nets
return nc
else:
return f(nets, *args, **kwargs)
else:
@wraps(f)
def wrapped_f(nets, *args, **kwargs):
if isinstance(net, NetCollection):
for net in nets:
f(nets, *args, **kwargs)
else:
f(nets, *args, **kwargs)
return wrapped_f
[docs]@net_util
def surface_of_revolution(curve, axis=2, name="Surface of Revolution"):
"""
Generate a surface of revolution from a planar curve
Parameters
----------
curve : ddg.datastructures.nets.net.SmoothNet
Planar curve.
axis : 0, 1, 2
Axis of rotation.
name: str
Name of the surface.
Returns
-------
ddg.datastructures.nets.net.SmoothNet
Surface of revolution.
"""
def fct(u, v):
x = np.array([curve(u)[0] * np.cos(v), curve(u)[0] * np.sin(v)])
if axis == 2:
x = np.append(x, curve(u)[1])
else:
x = np.insert(x, axis, curve(u)[1])
return x
# Might need a change when Intervals are removed
intervals = [curve.domain.recover()]
intervals.append([0, 2 * np.pi, True])
return (fct, intervals)
[docs]@net_util
def cylinder(curve, axis=2, name="Cylinder"):
"""
Generate a cylinder from a planar base curve
Parameters
----------
curve : ddg.datastructures.nets.net.SmoothNet
Planar curve.
axis : 0, 1, 2
Axis of the cylinder
name: str
Name of the surface.
Returns
-------
ddg.datastructures.nets.net.SmoothNet
Cylinder.
"""
def fct(u, v):
if axis == 2:
x = np.append(curve(u), v)
else:
x = np.insert(curve(u), axis, v)
return x
# Might need a change when Intervals are removed
intervals = [curve.domain.recover()]
intervals.append([-np.inf, np.inf])
return (fct, intervals)
[docs]@net_util
def cone(curve, cone_vertex, name="Cone"):
"""
Generate a cone over a given cone.
Parameters
----------
curve : ddg.datastructures.nets.net.SmoothNet
Spatial curve.
cone_vertex : numpy.array
Coordinates of a point in space
name: str
Name of the surface.
Returns
-------
ddg.datastructures.nets.net.SmoothNet
Cone
Notes
-----
* The cone is generated by connecting all the pointes of the curve
with the cone vertex.
"""
def fct(u, v):
x = np.multiply(curve(u), v) + np.multiply(cone_vertex, 1 - v)
return x
intervals = [curve.domain.recover()]
intervals.append([-np.inf, np.inf])
return (fct, intervals)
[docs]def sample_interval(interval, sample, option="", atol=None, anchor=None):
"""
Sample an interval.
Parameters
----------
interval : tuple
tuple containing the boundaries of the interval
sample : float, int or list
list is containing a float in the first, and int in the second entry
Stepsize (or amount of samples if 't' option is set). If given a list
for compound sampling, the first entry is the stepsize, while the
second is the total amount of samples.
option : string (default="")
| Options for the sampling.
|
| The supported options are:
| 't'otal amount of samples given
| 's'ymmetric sampling
| 'p'eriodic
| 'c'ompound sampling, both stepsize and amount of samples given
|
| When option 'c' is set, the sampling will choose between the stepsize
| and total amount according to the type of interval. In the case of an
| unbounded interval, it will also cut off the interval similar to
| bound_domain.
|
| When option 's' is set, the sampling will be fit symmetrically inside
| the interval.
atol : float or None (default=None)
Tolerance for the sampling. Is not used when 't' is set or the
interval is unbounded.
This function uses the global tolerance defaults if `atol` or `rtol`
are set to None. See :py:mod:`ddg.nonexact` for details.
anchor : float or None (default=None)
Point inside the interval.
The sampling will be chosen such that the point is part of it.
Returns
-------
tuple
| containing:
| [0] the (modified) boundaries of the sampled interval
| [1] amount of samples
| [2] sampling function
Raises
------
NotImplementedError
if given options are not compatible with eachother
ValueError
if sample is not compatible with the given options
See Also
--------
bound_domain
"""
a, b = interval
# Check if interval is bounded
# TODO: Total amounts and anchor given
valid_options = set(["s", "t", "p", "c"])
if not set(option).issubset(valid_options):
raise NotImplementedError(
"Unkown options found: " + str(set(option).difference(valid_options))
)
if "t" in option and "s" in option:
raise NotImplementedError("Uncompatible options 's' and 't'.")
if "t" in option and "c" in option:
raise NotImplementedError("Uncompatible options 'c' and 't'.")
if "c" in option and "s" in option:
raise NotImplementedError("Uncompatible options 'c' and 's'.")
if "c" in option and len(sample) != 2:
raise ValueError("Compound sampling needs to be [stepsize,samples,c].")
if "c" not in option and isinstance(sample, Iterable) and len(sample) != 1:
raise ValueError("Only one value supported for non-compound sampling.")
if "c" not in option and isinstance(sample, Iterable):
sample = sample[0]
if is_bounded_interval(interval):
# Total amount of samples given
if "c" in option:
sample = sample[1]
option += "t"
if "t" in option and not anchor:
N = sample
if not isinstance(N, int):
raise ValueError("Total samples must be an integer.")
m, n = 0, N
if "p" in option:
N += 1
def tmp(x, a=a, b=b, m=m, N=N):
return (x - m) * (b - a) / (N - 1) + a
return ([a, b], n, tmp)
elif "t" in option:
msg = "Total samples is not compatible with anchor."
raise NotImplementedError(msg)
# Stepsize given
else:
# Anchor given
if anchor is not None:
if "s" in option:
msg = "Option 's' is not compatible with anchor."
raise NotImplementedError(msg)
else:
N, re = divmod(anchor - a, sample)
N = int(np.floor(N))
if not nonexact.isclose(sample, re, atol=atol):
a += re
tempsample = np.arange(a, b + sample, sample)
if tempsample[-1] > b:
tempsample = tempsample[:-1]
if nonexact.isclose(tempsample[-1], b, atol=atol):
tempsample[-1] = b
N = len(tempsample)
m, n = 0, N
if "p" in option and tempsample[-1] == b:
N += 1
def tmp(x, a=a, b=tempsample[-1], m=m, N=N):
return (x - m) * (b - a) / (N - 1) + a
return ([a, tempsample[-1]], n, tmp)
# Symmetric case
elif "s" in option:
N, re = divmod((b - a), sample)
N = int(np.floor(N))
if nonexact.isclose(re, sample, atol=atol):
N += 2
else:
N += 1
if not nonexact.isclose(re, 0, atol=atol, rtol=0.0):
b = a + re / 2 + (N - 1) * sample
a += re / 2
m, n = 0, N
def tmp(x, a=a, b=b, m=m, N=N):
return (x - m) * (b - a) / (N - 1) + a
return ([a, b], n, tmp)
# No options given
else:
tempsample = np.arange(a, b + sample, sample)
if tempsample[-1] > b:
tempsample = tempsample[:-1]
if nonexact.isclose(tempsample[-1], b, atol=atol):
tempsample[-1] = b
N = len(tempsample)
m, n = 0, N
def tmp(x, a=a, b=tempsample[-1], m=m, N=N):
return (x - m) * (b - a) / (N - 1) + a
return ([a, tempsample[-1]], n, tmp)
else:
total = None
if "c" in option:
total = sample[1]
sample = sample[0]
if "t" in option:
msg = "Option 't' is not compatible with unbounded intervals."
raise ValueError(msg)
if a == -np.inf and b == np.inf:
if total:
N, re = divmod(total, 2)
b = N * sample
a = (1 - N - re) * sample
def tmp(x, s=sample, a=1 - N - re):
return s * (x + a)
else:
def tmp(x, s=sample):
return s * x
elif a == -np.inf and b < np.inf:
if total:
a = b - sample * (total - 1)
def tmp(x, s=sample, a=b - total + 1):
return s * (x + a)
else:
def tmp(x, s=sample, b=b):
return s * x + b
elif b == np.inf:
def tmp(x, s=sample, a=a):
return s * x + a
if total:
b = a + sample * (total - 1)
else:
raise ValueError("Strange Interval given")
if not total:
total = np.inf
return ([a, b], total, tmp)
[docs]def coordinate_hypersurface(net, component, value=None):
"""Generate the (parameterized) level set of a coordinate-direction.
Parameters
----------
net : SmoothNet or DiscreteNet with rectangular domain
Source net.
component : numpy.ndarray
Coordinate to set to a constant value.
value : Value for the chosen component of the domain.
If no value is given the output net will depend on the value as a Parameter.
Returns
-------
ddg.datastructures.nets.net.Net
The hypersurface as a net.
Notes
-----
* The hypersurface is generated by setting one coordinates to a
constant given value.
"""
if value is None:
default_value = (net.domain[component][0] + net.domain[component][1]) / 2
intervals = (
net.domain.recover()[:component] + net.domain.recover()[component + 1 :]
)
domain = net.domain.__class__(intervals)
def fct(*args, value=0.0, **kwargs):
args = args[:component] + (value,) + args[component:]
return net(*args, **kwargs)
if isinstance(net.domain, SmoothRectangularDomain):
new_net = SmoothNet(fct, domain)
elif isinstance(net.domain, DiscreteRectangularDomain):
new_net = DiscreteNet(fct, domain)
else:
intervals = (
net.domain.recover()[:component] + net.domain.recover()[component + 1 :]
)
domain = net.domain.__class__(intervals)
def fct(*args, **kwargs):
args = args[:component] + (value,) + args[component:]
return net(*args, **kwargs)
if isinstance(net.domain, SmoothRectangularDomain):
new_net = SmoothNet(fct, domain)
elif isinstance(net.domain, DiscreteRectangularDomain):
new_net = DiscreteNet(fct, domain)
return new_net
[docs]def coordinate_line(net, point, direction):
"""Gives a curve on the net starting at a point in a given direction.
Parameters
----------
net : SmoothNet or DiscreteNet with rectangular domain
Given net.
point : Iterable of Real
Point in the domain from where to start the curve.
direction : numpy.array, or int
Direction in the domain as an array, or component as int.
Returns
-------
ddg.datastructures.nets.net.Net
A line on the net.
Notes
-----
* If direction is given as an int, the periodicity of the net in the given
direction is inherited by the curve.
"""
point = np.array(point, dtype=np.float64)
if isinstance(direction, int):
direction = np.eye(net.dimension)[direction]
direction = np.array(direction, dtype=np.float64)
periodic = False
if isinstance(direction, int):
periodic = direction in net.domain.periodicity
bounds = []
for v, (a, b), co in zip(direction, net.domain.intervals, point):
if v != 0:
if v < 0:
a, b = b, a
bounds.append([(a - co) / v, (b - co) / v])
bounds = np.transpose(bounds)
intervals = [[np.nanmax(bounds[0]), np.nanmin(bounds[1]), periodic]]
domain = net.domain.__class__(intervals)
fct = lambda s: net(*(point + (s * direction)))
if isinstance(net, SmoothNet):
coordinate_line = SmoothCurve(fct, domain)
elif isinstance(net, DiscreteNet):
coordinate_line = DiscreteCurve(fct, domain)
else:
raise ValueError("Only works for smooth or discrete nets.")
return coordinate_line
[docs]def coordinate_lines(net, direction, sampling, anchor=None, atol=None):
"""Get coordinate lines in a given coordinate direction.
Parameters
----------
net : ddg.datastructures.nets.net.SmoothNet
Given net (must be SmoothNet with rectangular domain).
direction : int
Direction of the coordinate lines (component number).
sampling : list, int or float
See sample_smooth_domain
anchor : list, tuple or NONE
Point inside the domain of the net.
The sampling will be chosen such that the point is part of it.
atol : float or None (default=None)
Is used as tolerance for symmetric sampling.
See: sample_interval
This function uses the global tolerance defaults if `atol` or `rtol`
are set to None. See :py:mod:`ddg.nonexact` for details.
Returns
-------
ddg.datastructures.nets.net.NetCollection
Coordinate lines as a NetCollection.
Notes
-----
* They will be stepsize apart in all other directions.
* A set anchor precedes symmetric sampling
"""
if not isinstance(direction, int):
raise ValueError("Direcion should be given by one int value.")
anchor = np.array(anchor)
dom = (
net.domain.recover()
) # list of intervals with periodictiy information, len(dom)= dimension
direction_value = dom[direction][
0
] # one point of the domain, we need it for embedding
small_domain = delete_direction(net.domain, direction)
if anchor is not None:
anchor = np.delete(anchor, direction)
sampled_domain = netsconv.sample_smooth_domain(
small_domain, sampling, anchor=anchor, atol=atol
)
sampled_domain = embed(sampled_domain, level=direction_value, component=direction)
nc = NetCollection()
for point in sampled_domain:
cl = coordinate_line(net, point, direction)
nc.add(cl)
return nc
[docs]def coordinate_grid(net, sampling, anchor=None, atol=None):
"""Get a grid of your net of given sampling (fineness in each direction)
Parameters
----------
net: ddg.datastructures.nets.net.SmoothNet
Given net (must be SmoothNet with RectangularDomain).
sampling : list, int or float
Determines how the domain is sampled. See sample_smooth_domain.
anchor : numpy.ndarray, or None
Point inside the domain of the net.
The sampling will be chosen such that the point is part of it.
atol : float or None (default=None)
Is used as tolerance for symmetric sampling.
See: sample_interval
This function uses the global tolerance defaults if `atol` or `rtol`
are set to None. See :py:mod:`ddg.nonexact` for details.
Returns
-------
list of NetCollections of coordinate lines in every direction
"""
nclist = []
dim = net.domain.dimension
sampling_list = netsconv.sampling_decomposer(sampling, atol, anchor, dim)
for i in range(dim):
# this if else stuff is here because the sampling decomposer works a
# little weird. sometimes the stepsize is in a list, sometimes not
if isinstance(sampling_list[0][i], Iterable):
coordlines = coordinate_lines(
net,
i,
[sampling_list[0][i][0], sampling_list[1][i]],
anchor=sampling_list[3],
atol=sampling_list[2],
) # creates coordlines of net in direction i
else:
coordlines = coordinate_lines(
net,
i,
[sampling_list[0][i], sampling_list[1][i]],
anchor=sampling_list[3],
atol=sampling_list[2],
) # creates coordlines of net in direction i
nclist.append(coordlines)
return nclist
[docs]def diagonal_lines(net, direction, sampling, anchor=None, atol=None):
"""Get coordinate lines in a given coordinate direction.
Parameters
----------
net : ddg.datastructures.nets.net.SmoothNet
Given 2-dim net (must be SmoothNet with rectangular domain). Only works
for bounded nets.
direction: np.array of dim 2
vector consisting of 1 and -1, determining the direction of the diagonals
sampling : list of length 2 of form [number, "option"] or float (stepsize)
See ddg.datastructures.nets.conversion.sample_smooth_domain
anchor : numpy.ndarray, or None
Point inside the domain of the net.
The sampling will be chosen such that the point is part of it.
atol : float or None (default=None)
Is used as tolerance for symmetric sampling.
See: sample_interval
This function uses the global tolerance defaults if `atol` or `rtol` are
set to None. See :py:mod:`ddg.nonexact` for details.
Returns
-------
ddg.datastructures.nets.net.NetCollection
Coordinate lines as a NetCollection.
Notes
-----
* They will be stepsize apart in all other directions.
* A set anchor precedes symmetric sampling
"""
nc = NetCollection()
dom = (
net.domain.recover()
) # list of intervals with periodictiy information, len(dom)= dimension
if len(dom) != 2:
raise ValueError("Only works for 2-dim Nets")
# decoding our domain as [x1,x2]x[y1,y2]
x1 = dom[0][0]
x2 = dom[0][1]
y1 = dom[1][0]
y2 = dom[1][1]
sample_left = False
sample_right = False
# 1st case: diagonals in direction of vector [1,1]. we need to sample the
# left side of the rectangle
if np.array_equal(direction, [1, 1]) or np.array_equal(direction, [-1, -1]):
sample_left = True
# 2nd case: diagonals in direction [1,-1]. we need to sample the right side
# of the rectangle
elif np.array_equal(direction, [1, -1]) or np.array_equal(direction, [-1, 1]):
sample_right = True
else:
raise ValueError("Direction must consist of 1 and -1 values.")
bottom_side = delete_direction(net.domain, 1)
if sample_left:
extended_bottom_side = modify_direction(bottom_side, 0, [x1 - (y2 - y1), x2])
if isinstance(sampling, Iterable):
if len(sampling) != 2:
raise ValueError(
"You can only specify the sampling of one direction of your domain."
" Only works for bounded domains."
)
if anchor is not None:
raise NotImplementedError(
"Options 's' and 't' not compatible with anchor"
)
sampled_extended_bottom = netsconv.sample_smooth_domain(
extended_bottom_side, sampling, atol=atol
)
# Now we need to add the coordinate lines at the points in the
# first interval(bottom side), and the intersection points of the
# diags with the second interval (left side)
for idx in sampled_extended_bottom.domain.traverser:
bottom_point = sampled_extended_bottom(idx[0])[0]
if bottom_point < x1:
# we want to find the points on the "left side" we want to
# send coord lines from
embedded_point = [x1, x1 - bottom_point + y1]
cl = coordinate_line(net, embedded_point, direction)
nc.add(cl)
elif bottom_point >= x1: # we just need to embed the bottom points
embedded_point = [bottom_point, y1]
cl = coordinate_line(net, embedded_point, direction)
nc.add(cl)
elif isinstance(sampling, float) or isinstance(sampling, int):
if isinstance(sampling, int):
sampling = float(sampling)
if anchor is None:
bottom_anchor = x1
if anchor is not None:
bottom_anchor = anchor[0] - (anchor[1] - y1)
sampled_extended_bottom = netsconv.sample_smooth_domain(
extended_bottom_side, sampling, anchor=[bottom_anchor], atol=atol
)
# CALCULATING sampled points of LEFT side as intersection of the
# diaglines in extended bottom side
for idx in sampled_extended_bottom.domain.traverser:
bottom_point = sampled_extended_bottom(idx[0])[0]
if bottom_point < x1: # THEN WE WANT TO SAMPLE THE LEFT SIDE
embedded_point = [x1, x1 - bottom_point + y1]
cl = coordinate_line(net, embedded_point, direction)
nc.add(cl)
elif bottom_point >= x1:
embedded_point = [bottom_point, y1]
cl = coordinate_line(net, embedded_point, direction)
nc.add(cl)
else:
raise NotImplementedError(
"Only works for sampling options 's' and 't' or given float as stepsize"
)
if sample_right:
extended_bottom_side = modify_direction(bottom_side, 0, [x1, (y2 - y1) + x2])
if isinstance(sampling, Iterable):
if len(sampling) != 2:
raise ValueError(
"You can only specify the sampling of one direction of your domain."
" Only works for bounded domains."
)
if anchor is not None:
raise NotImplementedError(
"Options 's' and 't' not compatible with anchor"
)
sampled_extended_bottom = netsconv.sample_smooth_domain(
extended_bottom_side, sampling, atol=atol
)
# Now we need to add the coordinate lines at the points in the
# first interval(bottom side), and the intersection points of the
# diags with the second interval (right side)
for idx in sampled_extended_bottom.domain.traverser:
bottom_point = sampled_extended_bottom(idx[0])[0]
if bottom_point > x2:
# we want to find the points on the "right side" we want to
# send coord lines from
embedded_point = [x2, bottom_point - x2 + y1]
cl = coordinate_line(net, embedded_point, direction)
nc.add(cl)
elif bottom_point <= x2: # we just need to embed the bottom points
embedded_point = [bottom_point, y1]
cl = coordinate_line(net, embedded_point, direction)
nc.add(cl)
elif isinstance(sampling, float) or isinstance(sampling, int):
if isinstance(sampling, int):
sampling = float(sampling)
if anchor is None:
bottom_anchor = x2
if anchor is not None:
bottom_anchor = anchor[0] + (anchor[1] - y1)
sampled_extended_bottom = netsconv.sample_smooth_domain(
extended_bottom_side, sampling, anchor=[bottom_anchor], atol=atol
)
# CALCULATING sampled points of LEFT side as intersection of the
# diaglines in extended bottom side
for idx in sampled_extended_bottom.domain.traverser:
bottom_point = sampled_extended_bottom(idx[0])[0]
if bottom_point > x2:
# we want to find the points on the "right side" we want to
# send coord lines from
embedded_point = [x2, bottom_point - x2 + y1]
cl = coordinate_line(net, embedded_point, direction)
nc.add(cl)
elif bottom_point <= x2: # we just need to embed the bottom points
embedded_point = [bottom_point, y1]
cl = coordinate_line(net, embedded_point, direction)
nc.add(cl)
else:
raise NotImplementedError(
"Only works for sampling options 's' and 't' or given float as stepsize"
)
return nc
[docs]def octahedral_grid(net, stepsize, anchor=None, atol=None):
"""Get diagonal lines in each diagonal direction of the coordinate planes
of the domain
Parameters
----------
net : ddg.datastructures.nets.net.SmoothNet
Given 3-dimensional net (must be SmoothNet with rectangular domain).
stepsize : int or float
Determines the length of the diagonal of the occuring squares.
anchor : numpy.ndarray of dimension 3, or None
Point inside the domain of the net.
The sampling will be chosen such that the point is part of it.
atol : float or None (default=None)
Is used as tolerance for symmetric sampling.
See: sample_interval
This function uses the global tolerance defaults if `atol` or `rtol`
are set to None. See :py:mod:`ddg.nonexact` for details.
Returns
-------
list of NetCollections
ddg.datastructures.nets.net.NetCollection
Diagonal lines as a NetCollection.
"""
nc_xy = NetCollection()
nc_yz = NetCollection()
nc_xz = NetCollection()
dom = (
net.domain.recover()
) # list of intervals with periodictiy information, len(dom)= dimension
if len(dom) != 3:
raise ValueError("Only works for 3-dim Nets")
if anchor is None:
anchor = np.array([dom[0][0], dom[1][0], dom[2][0]])
if anchor is not None:
# FIRST DIRECTION
xy_plane = delete_direction(net.domain, 2)
# we sample the z-axis in half the stepsize
z_axis = delete_direction(net.domain, 0)
z_axis = delete_direction(
z_axis, 0
) # we delete the first and the 2nd direction
sampled_z_axis = netsconv.sample_smooth_domain(
z_axis, stepsize / 2, anchor=(anchor[2],), atol=atol
)
c = 0
for point in sampled_z_axis:
def new_fct_2d(a, b, p=point[0], **kwargs):
# a point of the sampled z axis is a list of size one,
# therefore point[0]
return net.fct(a, b, p, **kwargs)
new_net = SmoothNet(new_fct_2d, xy_plane)
if c % 2 == 0:
anch_2d = np.array([anchor[0], anchor[1]])
else:
s = stepsize / 2
if anchor[0] - s >= dom[0][0]:
anch_2d = np.array([anchor[0] - s, anchor[1]])
else:
anch_2d = np.array([anchor[0] + s, anchor[1]])
for line in diagonal_lines(
new_net, [1, 1], stepsize, anchor=anch_2d, atol=atol
):
nc_xy.add(line)
for line in diagonal_lines(
new_net, [-1, 1], stepsize, anchor=anch_2d, atol=atol
):
nc_xy.add(line)
c += 1
# SECOND DIRECTION
xz_plane = delete_direction(net.domain, 1)
# sampling the y axis
y_axis = delete_direction(net.domain, 0)
y_axis = delete_direction(y_axis, 1)
sampled_y_axis = netsconv.sample_smooth_domain(
y_axis, stepsize / 2, anchor=(anchor[1],), atol=atol
)
c = 0
for point in sampled_y_axis:
def new_fct_2d(a, b, p=point[0], **kwargs):
return net.fct(a, p, b, **kwargs)
new_net = SmoothNet(new_fct_2d, xz_plane)
if c % 2 == 0:
anch_2d = np.array([anchor[0], anchor[2]])
else:
s = stepsize / 2
if anchor[0] - s >= dom[0][0]:
anch_2d = np.array([anchor[0] - s, anchor[2]])
else:
anch_2d = np.array([anchor[0] + s, anchor[2]])
for line in diagonal_lines(
new_net, [1, 1], stepsize, anchor=anch_2d, atol=atol
):
nc_xz.add(line)
for line in diagonal_lines(
new_net, [-1, 1], stepsize, anchor=anch_2d, atol=atol
):
nc_xz.add(line)
c += 1
# THIRD DIRECTION
yz_plane = delete_direction(net.domain, 0)
# sampling the x axis
x_axis = delete_direction(net.domain, 1)
x_axis = delete_direction(x_axis, 1) # deleting second and third dir
sampled_x_axis = netsconv.sample_smooth_domain(
x_axis, stepsize / 2, anchor=(anchor[0],), atol=atol
)
c = 0
for point in sampled_x_axis:
def new_fct_2d(a, b, p=point[0], **kwargs):
return net.fct(p, a, b, **kwargs)
new_net = SmoothNet(new_fct_2d, yz_plane)
if c % 2 == 0:
anch_2d = np.array([anchor[1], anchor[2]])
else:
s = stepsize / 2
if anchor[0] - s >= dom[0][0]:
anch_2d = np.array([anchor[1] - s, anchor[2]])
else:
anch_2d = np.array([anchor[1] + s, anchor[2]])
for line in diagonal_lines(
new_net, [1, 1], stepsize, anchor=anch_2d, atol=atol
):
nc_yz.add(line)
for line in diagonal_lines(
new_net, [-1, 1], stepsize, anchor=anch_2d, atol=atol
):
nc_yz.add(line)
c += 1
return [nc_yz, nc_xz, nc_xy]
[docs]@domain_util
def shrink_domain(domain, amount):
"""
Shrinks domain by a given amount from the left and right.
Parameters
----------
domain : ddg.datastructures.nets.domain.DiscreteRectangularDomain or SmoothRectangularDomain
Domain to shrink, or a net to take the domain from.
amount : list or float
| Specifies how much is taken away in each direction from the domain:
| * float - value to take away from left and right in every direction
| * list - every entry in the list corresponds to one direction,
| and may be None, int, or a tuple of int or None.
Returns
-------
ddg.datastructures.nets.domain.Domain
Shrunken domain, if the domain was given.
ddg.datastructures.nets.net.Net
The input discrete net containing the new shrunken domain,
if the net was given.
Notes
-----
* Only bounded directions will be shrunken.
* Periodicity is inherited only for the directions that are not shrunken.
""" # noqa: E501
if not isinstance(amount, Iterable):
amount = [amount for _ in range(domain.dimension)]
shrunken_intervals = []
for i, (s, interval) in enumerate(zip(amount, domain.recover())):
if s is None or not is_bounded_interval(interval):
shrunken_intervals.append(interval)
else:
if len(np.shape(s)) == 0:
s = [s, s]
a, b, p = interval
s0, s1 = s
if s0 is not None:
a = a + s0
if s1 is not None:
b = b - s1
shrunken_intervals.append([a, b, p])
return shrunken_intervals
[docs]@domain_util
def bound_domain(domain, bounding):
"""
Parameters
----------
domain : ddg.datastructures.nets.domain.Smooth/DiscreteRectangularDomain or ddg.datastructures.nets.net.Smooth/DiscreteNet
Given (unbounded) rectangular domain to bound, or net
whose domain is to be bound
bounding : int/float
| Bounding range.
| If an interval (a,b) in the domain is unbounded,
| we bound it in the smooth case by
| (1) [a, a+bounding], if a > -inf
| (2) [b-bounding, b], if b < inf, or
| (3) [-bounding/2, bounding/2] if a == -inf and b == inf.
|
| In the discrete case, bounding is the total amount of points chosen in
| an unbounded direction.
| (1) [a, a+bounding-1], if a > -inf
| (2) [b-bounding+1, b], if b < inf, or
| (3) [1-bounding/2, bounding/2] if a == -inf and b == inf.
Returns
-------
ddg.datastructures.nets.domain.Smooth/DiscreteRectangularDomain or ddg.datastructures.nets.net.Smooth/DiscreteNet
Will be the same type as domain, e.g. a discrete/smooth
rectangular domain, if domain was discrete/smooth rectangular domain.
""" # noqa: E501
# Empty domains get handled by the decorator
if domain.bounded:
# domain.recover() doesn't work because of SmoothInterval and
# DiscreteInterval
return [d + [i in domain.periodicity] for i, d in enumerate(domain.intervals)]
discrete = 0
if isinstance(domain, DiscreteDomain):
discrete = 1
intervals = []
for idx, el in enumerate(domain):
if is_bounded_interval(el):
intervals.append(el + [idx in domain.periodicity])
else:
a, b = el
if a > -np.inf:
intervals.append([a, a + bounding - discrete])
elif b < np.inf:
intervals.append([b - bounding + discrete, b])
else:
n, re = divmod(bounding, 2)
intervals.append([discrete - re - n, n])
return intervals
[docs]@domain_util
def create_subdomain(domain, subdomain):
"""
Create subdomain of a given discrete domain.
WARNING:
To assure that the output is always a subdomain
of the given domain, subdomain will be intersected
with the domain.
Parameters
----------
domain : ddg.datastructures.nets.domain.DiscreteRectangularDomain
Given domain to reduce to a subdomain, or net.
subdomain : list, or DiscreteDomain
| Each entry in the list corresponds to a direction:
| * int - number of values to take
| * tuple - subdomain in given direction
| If a DiscreteDomain is given, it will just be returned as is.
Returns
-------
ddg.datastrucuters.nets.domain.DiscreteRectangularDomain
Subdomain.
ddg.datastructures.nets.net.Net
A discrete net containing the new shrunken domain,
if the net was given.
"""
intervals = []
if subdomain is None:
subdomain = [None] * domain.dimension
for idx, (el, interval) in enumerate(zip(subdomain, domain)):
if el is None or el == interval:
intervals.append(interval + [idx in domain.periodicity])
else:
ela, elb = el
ia, ib = interval
intervals.append([max(ela, ia), min(elb, ib)])
return intervals
[docs]@net_util
def continue_by_reflection(net, direction, reflection, invertdirection=False):
"""
Continue a net by reflection.
Given a net, direction and reflection, the function returns a net with an
expanded domain in direction direction. Points outside the domain of the
orignal domain will be mapped to a point corresponding to it inside of the
domain before the new net evaluates first the function of the original net
before applying the reflection.
Examples
--------
>>> import numpy as np
>>> from ddg.datastructures.nets.net import SmoothNet
>>> from ddg.datastructures.nets.utils import continue_by_reflection
>>> net = SmoothNet(lambda *a: np.array(a), [[0, 2]] * 2)
>>> reflectionnet = continue_by_reflection(net, 0, 0)
>>> reflectionnet(3, 0)
array([-1, 0])
Parameters
----------
net : ddg.datastructures.nets.net.SmoothNet or DiscreteNet
Net to continue by reflection
direction : int
Direction in which the domain will be expanded. Points outside the
original domain will be mapped to the point corresponding to it and the
reflection will be applied after the function was evalued
reflection : int, Iterable, function
| If given an integer or an Iterable containing integers, reflection
| determines the coordinates changing sign after the original net is
| evalued at the point corresponding to the one given outside the
| original domain, e.g. for a point p outside of the original domain
| corresponding to p~ and reflect the function given by the int/Iterable
|
| newNet(p) == reflect(net(p~))
|
| If a function is given, it will be used instead.
|
| WARNING:
| In case a function is given, it will not be checked if it can
| accept the output of the net or not.
invertdirection : bool, optional
If set to False, the domain will be expanded in positive direction
If set to True, the domain will be expanded in negative direction
Returns
-------
ddg.datastructures.nets.net.Net
Same type as net
Raises
------
IndexError
If direction is not in range of domain dimension
TypeError
If reflection is neither an int, Iterable or function
"""
domain = net.domain
if direction >= domain.dimension:
raise IndexError("Dimension of net domain is < {}.".format(direction + 1))
if not callable(reflection):
if isinstance(reflection, Iterable) and all(
[isinstance(x, int) for x in reflection]
):
def rflct(*args, idx=reflection):
nargs = np.array(args)
reflect = np.array([-1 if i in idx else 1 for i in range(len(nargs))])
return nargs * reflect
elif isinstance(reflection, int):
def rflct(*args, idx=reflection):
nargs = np.array(args)
reflect = np.array([-1 if i == idx else 1 for i in range(len(nargs))])
return nargs * reflect
else:
msg = "Does not support {} as reflection.".format(type(reflection))
raise TypeError(msg)
reflection = rflct
if not invertdirection:
intervals = []
for idx, el in enumerate(domain):
if idx != direction:
intervals.append(el + [idx in domain.periodicity])
else:
a, b = el
intervals.append([a, 2 * b - a])
def fct(*args, **kwargs):
bound = domain.intervals[direction][1]
if args in domain:
return net(*args, **kwargs)
else:
nargs = np.array(args)
nargs[direction] = 2 * bound - nargs[direction]
return reflection(*net(*nargs))
else:
intervals = []
for idx, el in enumerate(domain):
if idx != direction:
intervals.append(el + [idx in domain.periodicity])
else:
a, b = el
intervals.append([2 * a - b, b])
def fct(*args, **kwargs):
bound = domain.intervals[direction][0]
if args in domain:
return net(*args, **kwargs)
else:
nargs = np.array(args)
nargs[direction] = 2 * bound - nargs[direction]
return reflection(*net(*nargs))
return (fct, intervals)
[docs]@net_util
def concatenate(net, fct):
"""
Create a net which's function is the conatenation of the function of the
original net and fct and has the same domain as the original.
Parameters
----------
net : ddg.datastructures.nets.Net or subclass of it
Net with compatible function for fct
fct : function
Function to concatenate with the net function
Returns
-------
ddg.datastructures.nets.Net
Exact type is the same as given net
"""
def nfct(*args, **kwargs):
return fct(net(*args))
intervals = net.domain.intervals
return (nfct, intervals)
[docs]def evaluate(net):
"""
Evaluate a discrete net on its entire domain.
This should NOT replace net[net.domain]. Instead use this function if the
function of a net has effects outside it, i.e. when a net creates blender
object and link them to the scene on its own.
Parameters
----------
net : ddg.datastructures.nets.net.DiscreteNet or subclass of it
Returns
-------
0
If successfull
Raises
------
TypeError
if net was not a (subclass of) ddg.datastructures.nets.net.DiscreteNet
"""
if not isinstance(net, DiscreteNet):
raise TypeError("Only works for discrete nets.")
return net.evaluate()
[docs]@domain_util
def delete_direction(domain, direction):
"""
Delete one direction of a domain.
Parameters
----------
domain: SmoothRectangularDomain or Discrete RectangularDomain
may be parametrized
direction: int
Returns
-------
domain
of same class of one dimension less
"""
intdict = {}
idx = 0
for interval in domain: # creating dictionary -> index:interval
intdict[idx] = interval
idx += 1
del intdict[direction]
less_intervals = list(intdict.values()) # conv dictionary of intervals to list
return less_intervals
[docs]@domain_util
def modify_direction(domain, direction, new_interval):
intdict = {}
idx = 0
for interval in domain: # creating dictionary -> index:interval
intdict[idx] = interval
idx += 1
intdict[direction] = new_interval
new_intervals = list(intdict.values()) # conv dictionary of intervals to list
return new_intervals
[docs]@applicable_to_netcollection
def cut_bounding_box(curve, bounding_box=(np.inf, np.inf, np.inf)):
"""Cuts the curve within the given bounding box.
* only works for curves with infinite domain
* works best for curves with small sample size (no additional points on the
boundary of the bounding_box are set)
Parameters
----------
curve : ddg.DiscreteCurve
curve to be cut
bounding_box : iterable of three floats
distances of planes parallel to the coordinate planes
For each given entry the curve is cut in positive and negative direction.
Raises
------
ValueError
if the curves domain is infinite
Returns
-------
ddg.NetCollection
NetCollection of the remaining segments of the curve
"""
bounding_box = np.array([np.abs(x) for x in bounding_box])
if not is_bounded_interval(curve.domain.interval):
raise ValueError(
"The curves domain has to be bounded, the given curves domain is "
"unbounded:",
curve.domain.interval,
)
interval_collection = []
out = True
for i in range(curve.domain.interval[0], curve.domain.interval[1] + 1):
vec = np.array(curve.fct(i))
if (
not (np.array(vec) < bounding_box).all()
or not (np.array(vec) > -bounding_box).all()
):
out = True
else:
if out:
interval_collection.append([i, i])
out = False
else:
interval_collection[-1][-1] = i
curve_segments = []
if curve.domain.periodicity:
if len(interval_collection) >= 2:
I1 = interval_collection.pop(0)
I2 = interval_collection.pop()
if I1[0] == curve.domain.interval[0] and I2[1] == curve.domain.interval[1]:
I2[1] += I1[1] + 1
def periodicity_fct(x, b):
if x <= b:
return curve.fct(x)
else:
return curve.fct(x - b - 1)
curve_segments.append(ddg.DiscreteCurve(periodicity_fct, I2))
for interval in interval_collection:
curve_segments.append(ddg.DiscreteCurve(curve.fct, interval))
if len(curve_segments) == 1:
net = curve_segments[0]
else:
net = NetCollection(curve_segments)
return net