ISL -> isl
[linpy.git] / linpy / geometry.py
1 # Copyright 2014 MINES ParisTech
2 #
3 # This file is part of LinPy.
4 #
5 # LinPy is free software: you can redistribute it and/or modify
6 # it under the terms of the GNU General Public License as published by
7 # the Free Software Foundation, either version 3 of the License, or
8 # (at your option) any later version.
9 #
10 # LinPy is distributed in the hope that it will be useful,
11 # but WITHOUT ANY WARRANTY; without even the implied warranty of
12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 # GNU General Public License for more details.
14 #
15 # You should have received a copy of the GNU General Public License
16 # along with LinPy. If not, see <http://www.gnu.org/licenses/>.
17
18 import math
19 import numbers
20 import operator
21
22 from abc import ABC, abstractproperty, abstractmethod
23 from collections import OrderedDict, Mapping
24
25 from .linexprs import Symbol
26
27
28 __all__ = [
29 'GeometricObject',
30 'Point',
31 'Vector',
32 ]
33
34
35 class GeometricObject(ABC):
36 """
37 GeometricObject is an abstract class to represent objects with a
38 geometric representation in space. Subclasses of GeometricObject are
39 Polyhedron, Domain and Point.
40 """
41
42 @abstractproperty
43 def symbols(self):
44 """
45 The tuple of symbols present in the object expression, sorted according
46 to Symbol.sortkey().
47 """
48 pass
49
50 @property
51 def dimension(self):
52 """
53 The dimension of the object, i.e. the number of symbols present in it.
54 """
55 return len(self.symbols)
56
57 @abstractmethod
58 def aspolyhedron(self):
59 """
60 Return a Polyhedron object that approximates the geometric object.
61 """
62 pass
63
64 def asdomain(self):
65 """
66 Return a Domain object that approximates the geometric object.
67 """
68 return self.aspolyhedron()
69
70
71 class Coordinates:
72 """
73 This class represents coordinate systems.
74 """
75
76 __slots__ = (
77 '_coordinates',
78 )
79
80 def __new__(cls, coordinates):
81 """
82 Create a coordinate system from a dictionary or a sequence that maps the
83 symbols to their coordinates. Coordinates must be rational numbers.
84 """
85 if isinstance(coordinates, Mapping):
86 coordinates = coordinates.items()
87 self = object().__new__(cls)
88 self._coordinates = []
89 for symbol, coordinate in coordinates:
90 if not isinstance(symbol, Symbol):
91 raise TypeError('symbols must be Symbol instances')
92 if not isinstance(coordinate, numbers.Real):
93 raise TypeError('coordinates must be real numbers')
94 self._coordinates.append((symbol, coordinate))
95 self._coordinates.sort(key=lambda item: item[0].sortkey())
96 self._coordinates = OrderedDict(self._coordinates)
97 return self
98
99 @property
100 def symbols(self):
101 """
102 The tuple of symbols present in the coordinate system, sorted according
103 to Symbol.sortkey().
104 """
105 return tuple(self._coordinates)
106
107 @property
108 def dimension(self):
109 """
110 The dimension of the coordinate system, i.e. the number of symbols
111 present in it.
112 """
113 return len(self.symbols)
114
115 def coordinate(self, symbol):
116 """
117 Return the coordinate value of the given symbol. Raise KeyError if the
118 symbol is not involved in the coordinate system.
119 """
120 if not isinstance(symbol, Symbol):
121 raise TypeError('symbol must be a Symbol instance')
122 return self._coordinates[symbol]
123
124 __getitem__ = coordinate
125
126 def coordinates(self):
127 """
128 Iterate over the pairs (symbol, value) of coordinates in the coordinate
129 system.
130 """
131 yield from self._coordinates.items()
132
133 def values(self):
134 """
135 Iterate over the coordinate values in the coordinate system.
136 """
137 yield from self._coordinates.values()
138
139 def __bool__(self):
140 """
141 Return True if not all coordinates are 0.
142 """
143 return any(self._coordinates.values())
144
145 def __hash__(self):
146 return hash(tuple(self.coordinates()))
147
148 def __repr__(self):
149 string = ', '.join(['{!r}: {!r}'.format(symbol, coordinate)
150 for symbol, coordinate in self.coordinates()])
151 return '{}({{{}}})'.format(self.__class__.__name__, string)
152
153 def _map(self, func):
154 for symbol, coordinate in self.coordinates():
155 yield symbol, func(coordinate)
156
157 def _iter2(self, other):
158 if self.symbols != other.symbols:
159 raise ValueError('arguments must belong to the same space')
160 coordinates1 = self._coordinates.values()
161 coordinates2 = other._coordinates.values()
162 yield from zip(self.symbols, coordinates1, coordinates2)
163
164 def _map2(self, other, func):
165 for symbol, coordinate1, coordinate2 in self._iter2(other):
166 yield symbol, func(coordinate1, coordinate2)
167
168
169 class Point(Coordinates, GeometricObject):
170 """
171 This class represents points in space.
172
173 Point instances are hashable and should be treated as immutable.
174 """
175
176 def isorigin(self):
177 """
178 Return True if all coordinates are 0.
179 """
180 return not bool(self)
181
182 def __hash__(self):
183 return super().__hash__()
184
185 def __add__(self, other):
186 """
187 Translate the point by a Vector object and return the resulting point.
188 """
189 if isinstance(other, Vector):
190 coordinates = self._map2(other, operator.add)
191 return Point(coordinates)
192 return NotImplemented
193
194 def __sub__(self, other):
195 """
196 If other is a point, substract it from self and return the resulting
197 vector. If other is a vector, translate the point by the opposite vector
198 and returns the resulting point.
199 """
200 coordinates = []
201 if isinstance(other, Point):
202 coordinates = self._map2(other, operator.sub)
203 return Vector(coordinates)
204 elif isinstance(other, Vector):
205 coordinates = self._map2(other, operator.sub)
206 return Point(coordinates)
207 return NotImplemented
208
209 def __eq__(self, other):
210 """
211 Test whether two points are equal.
212 """
213 if isinstance(other, Point):
214 return self._coordinates == other._coordinates
215 return NotImplemented
216
217 def aspolyhedron(self):
218 from .polyhedra import Polyhedron
219 equalities = []
220 for symbol, coordinate in self.coordinates():
221 equalities.append(symbol - coordinate)
222 return Polyhedron(equalities)
223
224
225 class Vector(Coordinates):
226 """
227 This class represents vectors in space.
228
229 Vector instances are hashable and should be treated as immutable.
230 """
231
232 def __new__(cls, initial, terminal=None):
233 """
234 Create a vector from a dictionary or a sequence that maps the symbols to
235 their coordinates, or as the displacement between two points.
236 """
237 if not isinstance(initial, Point):
238 initial = Point(initial)
239 if terminal is None:
240 coordinates = initial._coordinates
241 else:
242 if not isinstance(terminal, Point):
243 terminal = Point(terminal)
244 coordinates = terminal._map2(initial, operator.sub)
245 return super().__new__(cls, coordinates)
246
247 def isnull(self):
248 """
249 Return True if all coordinates are 0.
250 """
251 return not bool(self)
252
253 def __hash__(self):
254 return super().__hash__()
255
256 def __add__(self, other):
257 """
258 If other is a point, translate it with the vector self and return the
259 resulting point. If other is a vector, return the vector self + other.
260 """
261 if isinstance(other, (Point, Vector)):
262 coordinates = self._map2(other, operator.add)
263 return other.__class__(coordinates)
264 return NotImplemented
265
266 def __sub__(self, other):
267 """
268 If other is a point, substract it from the vector self and return the
269 resulting point. If other is a vector, return the vector self - other.
270 """
271 if isinstance(other, (Point, Vector)):
272 coordinates = self._map2(other, operator.sub)
273 return other.__class__(coordinates)
274 return NotImplemented
275
276 def __neg__(self):
277 """
278 Return the vector -self.
279 """
280 coordinates = self._map(operator.neg)
281 return Vector(coordinates)
282
283 def __mul__(self, other):
284 """
285 Multiplies a Vector by a scalar value.
286 """
287 if isinstance(other, numbers.Real):
288 coordinates = self._map(lambda coordinate: other * coordinate)
289 return Vector(coordinates)
290 return NotImplemented
291
292 __rmul__ = __mul__
293
294 def __truediv__(self, other):
295 """
296 Divide the vector by the specified scalar and returns the result as a
297 vector.
298 """
299 if isinstance(other, numbers.Real):
300 coordinates = self._map(lambda coordinate: coordinate / other)
301 return Vector(coordinates)
302 return NotImplemented
303
304 def __eq__(self, other):
305 """
306 Test whether two vectors are equal.
307 """
308 if isinstance(other, Vector):
309 return self._coordinates == other._coordinates
310 return NotImplemented
311
312 def angle(self, other):
313 """
314 Retrieve the angle required to rotate the vector into the vector passed
315 in argument. The result is an angle in radians, ranging between -pi and
316 pi.
317 """
318 if not isinstance(other, Vector):
319 raise TypeError('argument must be a Vector instance')
320 cosinus = self.dot(other) / (self.norm()*other.norm())
321 return math.acos(cosinus)
322
323 def cross(self, other):
324 """
325 Compute the cross product of two 3D vectors. If either one of the
326 vectors is not three-dimensional, a ValueError exception is raised.
327 """
328 if not isinstance(other, Vector):
329 raise TypeError('other must be a Vector instance')
330 if self.dimension != 3 or other.dimension != 3:
331 raise ValueError('arguments must be three-dimensional vectors')
332 if self.symbols != other.symbols:
333 raise ValueError('arguments must belong to the same space')
334 x, y, z = self.symbols
335 coordinates = []
336 coordinates.append((x, self[y]*other[z] - self[z]*other[y]))
337 coordinates.append((y, self[z]*other[x] - self[x]*other[z]))
338 coordinates.append((z, self[x]*other[y] - self[y]*other[x]))
339 return Vector(coordinates)
340
341 def dot(self, other):
342 """
343 Compute the dot product of two vectors.
344 """
345 if not isinstance(other, Vector):
346 raise TypeError('argument must be a Vector instance')
347 result = 0
348 for symbol, coordinate1, coordinate2 in self._iter2(other):
349 result += coordinate1 * coordinate2
350 return result
351
352 def __hash__(self):
353 return super().__hash__()
354
355 def norm(self):
356 """
357 Return the norm of the vector.
358 """
359 return math.sqrt(self.norm2())
360
361 def norm2(self):
362 """
363 Return the squared norm of the vector.
364 """
365 result = 0
366 for coordinate in self._coordinates.values():
367 result += coordinate ** 2
368 return result
369
370 def asunit(self):
371 """
372 Return the normalized vector, i.e. the vector of same direction but with
373 norm 1.
374 """
375 return self / self.norm()