New example: Menger sponge, to be plotted
[linpy.git] / pypol / geometry.py
1 import math
2 import numbers
3 import operator
4
5 from abc import ABC, abstractproperty, abstractmethod
6 from collections import OrderedDict, Mapping
7
8 from .linexprs import Symbol
9
10
11 __all__ = [
12 'GeometricObject',
13 'Point',
14 'Vector',
15 ]
16
17
18 class GeometricObject(ABC):
19
20 @abstractproperty
21 def symbols(self):
22 pass
23
24 @property
25 def dimension(self):
26 return len(self.symbols)
27
28 @abstractmethod
29 def aspolyhedron(self):
30 pass
31
32 def asdomain(self):
33 return self.aspolyhedron()
34
35
36 class Coordinates:
37
38 __slots__ = (
39 '_coordinates',
40 )
41
42 def __new__(cls, coordinates):
43 if isinstance(coordinates, Mapping):
44 coordinates = coordinates.items()
45 self = object().__new__(cls)
46 self._coordinates = OrderedDict()
47 for symbol, coordinate in sorted(coordinates,
48 key=lambda item: item[0].sortkey()):
49 if not isinstance(symbol, Symbol):
50 raise TypeError('symbols must be Symbol instances')
51 if not isinstance(coordinate, numbers.Real):
52 raise TypeError('coordinates must be real numbers')
53 self._coordinates[symbol] = coordinate
54 return self
55
56 @property
57 def symbols(self):
58 return tuple(self._coordinates)
59
60 @property
61 def dimension(self):
62 return len(self.symbols)
63
64 def coordinates(self):
65 yield from self._coordinates.items()
66
67 def coordinate(self, symbol):
68 if not isinstance(symbol, Symbol):
69 raise TypeError('symbol must be a Symbol instance')
70 return self._coordinates[symbol]
71
72 __getitem__ = coordinate
73
74 def __bool__(self):
75 return any(self._coordinates.values())
76
77 def __hash__(self):
78 return hash(tuple(self.coordinates()))
79
80 def __repr__(self):
81 string = ', '.join(['{!r}: {!r}'.format(symbol, coordinate)
82 for symbol, coordinate in self.coordinates()])
83 return '{}({{{}}})'.format(self.__class__.__name__, string)
84
85 def _map(self, func):
86 for symbol, coordinate in self.coordinates():
87 yield symbol, func(coordinate)
88
89 def _iter2(self, other):
90 if self.symbols != other.symbols:
91 raise ValueError('arguments must belong to the same space')
92 coordinates1 = self._coordinates.values()
93 coordinates2 = other._coordinates.values()
94 yield from zip(self.symbols, coordinates1, coordinates2)
95
96 def _map2(self, other, func):
97 for symbol, coordinate1, coordinate2 in self._iter2(other):
98 yield symbol, func(coordinate1, coordinate2)
99
100
101 class Point(Coordinates, GeometricObject):
102 """
103 This class represents points in space.
104 """
105
106 def isorigin(self):
107 return not bool(self)
108
109 def __hash__(self):
110 return super().__hash__()
111
112 def __add__(self, other):
113 if not isinstance(other, Vector):
114 return NotImplemented
115 coordinates = self._map2(other, operator.add)
116 return Point(coordinates)
117
118 def __sub__(self, other):
119 coordinates = []
120 if isinstance(other, Point):
121 coordinates = self._map2(other, operator.sub)
122 return Vector(coordinates)
123 elif isinstance(other, Vector):
124 coordinates = self._map2(other, operator.sub)
125 return Point(coordinates)
126 else:
127 return NotImplemented
128
129 def __eq__(self, other):
130 return isinstance(other, Point) and \
131 self._coordinates == other._coordinates
132
133 def aspolyhedron(self):
134 from .polyhedra import Polyhedron
135 equalities = []
136 for symbol, coordinate in self.coordinates():
137 equalities.append(symbol - coordinate)
138 return Polyhedron(equalities)
139
140
141 class Vector(Coordinates):
142 """
143 This class represents displacements in space.
144 """
145
146 def __new__(cls, initial, terminal=None):
147 if not isinstance(initial, Point):
148 initial = Point(initial)
149 if terminal is None:
150 coordinates = initial._coordinates
151 else:
152 if not isinstance(terminal, Point):
153 terminal = Point(terminal)
154 coordinates = terminal._map2(initial, operator.sub)
155 return super().__new__(cls, coordinates)
156
157 def isnull(self):
158 return not bool(self)
159
160 def __hash__(self):
161 return super().__hash__()
162
163 def __add__(self, other):
164 if isinstance(other, (Point, Vector)):
165 coordinates = self._map2(other, operator.add)
166 return other.__class__(coordinates)
167 return NotImplemented
168
169 def angle(self, other):
170 """
171 Retrieve the angle required to rotate the vector into the vector passed
172 in argument. The result is an angle in radians, ranging between -pi and
173 pi.
174 """
175 if not isinstance(other, Vector):
176 raise TypeError('argument must be a Vector instance')
177 cosinus = self.dot(other) / (self.norm()*other.norm())
178 return math.acos(cosinus)
179
180 def cross(self, other):
181 """
182 Calculate the cross product of two Vector3D structures.
183 """
184 if not isinstance(other, Vector):
185 raise TypeError('other must be a Vector instance')
186 if self.dimension != 3 or other.dimension != 3:
187 raise ValueError('arguments must be three-dimensional vectors')
188 if self.symbols != other.symbols:
189 raise ValueError('arguments must belong to the same space')
190 x, y, z = self.symbols
191 coordinates = []
192 coordinates.append((x, self[y]*other[z] - self[z]*other[y]))
193 coordinates.append((y, self[z]*other[x] - self[x]*other[z]))
194 coordinates.append((z, self[x]*other[y] - self[y]*other[x]))
195 return Vector(coordinates)
196
197 def __truediv__(self, other):
198 """
199 Divide the vector by the specified scalar and returns the result as a
200 vector.
201 """
202 if not isinstance(other, numbers.Real):
203 return NotImplemented
204 coordinates = self._map(lambda coordinate: coordinate / other)
205 return Vector(coordinates)
206
207 def dot(self, other):
208 """
209 Calculate the dot product of two vectors.
210 """
211 if not isinstance(other, Vector):
212 raise TypeError('argument must be a Vector instance')
213 result = 0
214 for symbol, coordinate1, coordinate2 in self._iter2(other):
215 result += coordinate1 * coordinate2
216 return result
217
218 def __eq__(self, other):
219 return isinstance(other, Vector) and \
220 self._coordinates == other._coordinates
221
222 def __hash__(self):
223 return hash(tuple(self.coordinates()))
224
225 def __mul__(self, other):
226 if not isinstance(other, numbers.Real):
227 return NotImplemented
228 coordinates = self._map(lambda coordinate: other * coordinate)
229 return Vector(coordinates)
230
231 __rmul__ = __mul__
232
233 def __neg__(self):
234 coordinates = self._map(operator.neg)
235 return Vector(coordinates)
236
237 def norm(self):
238 return math.sqrt(self.norm2())
239
240 def norm2(self):
241 result = 0
242 for coordinate in self._coordinates.values():
243 result += coordinate ** 2
244 return result
245
246 def asunit(self):
247 return self / self.norm()
248
249 def __sub__(self, other):
250 if isinstance(other, (Point, Vector)):
251 coordinates = self._map2(other, operator.sub)
252 return other.__class__(coordinates)
253 return NotImplemented