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()