This chapter is taken from the book A Primer on Scientific Programming with Python by H. P. Langtangen, 5th edition, Springer, 2016.
This section explains how to implement two-dimensional vectors in Python such that these vectors act as objects we can add, subtract, form inner products with, and do other mathematical operations on. To understand the forthcoming material, it is necessary to have digested the section Special methods, in particular the sections Adding objects and Arithmetic operations and other special methods.
Vectors in the plane are described by a pair of real numbers, \( (a, b) \). There are mathematical rules for adding and subtracting vectors, multiplying two vectors (the inner or dot or scalar product), the length of a vector, and multiplication by a scalar: $$ \begin{align} (a,b) + (c,d) &= (a+c, b+d), \tag{4}\\ (a,b) - (c,d) &= (a-c, b-d), \tag{5}\\ (a,b)\cdot(c,d) &= ac + bd, \tag{6}\\ ||(a,b)|| &= \sqrt{(a,b)\cdot(a,b)} \thinspace . \tag{7} \end{align} $$ Moreover, two vectors \( (a,b) \) and \( (c,d) \) are equal if \( a=c \) and \( b=d \).
We may create a class for plane vectors where the above mathematical
operations are implemented by special methods. The class must contain
two data attributes, one for each component of the vector, called x and
y
below. We include special methods for addition, subtraction, the scalar
product (multiplication), the absolute value (length), comparison of
two vectors (== and !=), as well as a method for printing out
a vector.
class Vec2D(object):
def __init__(self, x, y):
self.x = x
self.y = y
def __add__(self, other):
return Vec2D(self.x + other.x, self.y + other.y)
def __sub__(self, other):
return Vec2D(self.x - other.x, self.y - other.y)
def __mul__(self, other):
return self.x*other.x + self.y*other.y
def __abs__(self):
return math.sqrt(self.x**2 + self.y**2)
def __eq__(self, other):
return self.x == other.x and self.y == other.y
def __str__(self):
return '(%g, %g)' % (self.x, self.y)
def __ne__(self, other):
return not self.__eq__(other) # reuse __eq__
The __add__,
__sub__,
__mul__,
__abs__, and
__eq__ methods
should be quite straightforward to understand from the previous
mathematical definitions of these operations.
The last method deserves a comment: here we simply reuse the
equality operator __eq__, but precede it with a not.
We could also have implemented this method as
def __ne__(self, other):
return self.x != other.x or self.y != other.y
Nevertheless, this implementation requires us to write more, and it has
the danger of introducing an error in the logics of the boolean
expressions. A more reliable approach, when we know that the
__eq__ method works, is to reuse this method and
observe that a != b means not (a == b).
A word of warning is in place regarding our implementation of the
equality operator (== via __eq__).
We test for equality of each component, which is correct from a
mathematical point of view. However, each vector component is a
floating-point number that may be subject to rounding errors both
in the representation on the computer and from previous (inexact)
floating-point calculations. Two mathematically equal components may
be different in their inexact representations on the computer.
The remedy for this problem is to avoid testing for equality, but instead
check that the difference between the components is sufficiently small.
The function numpy.allclose can be used for this purpose:
if a == b:
by
if numpy.allclose(a, b):
A more reliable equality operator can now be implemented:
class Vec2D(object):
...
def __eq__(self, other):
return numpy.allclose(self.x, other.x) and \
numpy.allclose(self.y, other.y)
As a rule of thumb, you should never apply the == test to two
float objects.
The special method
__len__ could be introduced as a synonym for
__abs__, i.e., for a Vec2D instance named v,
len(v) is the same as abs(v), because the absolute value of
a vector is mathematically the same as the length of the vector.
However, if we implement
def __len__(self):
# Reuse implementation of __abs__
return abs(self) # equiv. to self.__abs__()
we will run into trouble when we compute len(v) and the answer
is (as usual) a float. Python will then complain and tell us that
len(v) must return an int. Therefore,
__len__ cannot be used as a synonym for the length of the
vector in our application. On the other hand, we could let
len(v) mean the number of components in the vector:
def __len__(self):
return 2
This is not a very useful function, though, as we already know that all our
Vec2D vectors have just two components.
For generalizations of the class to vectors with \( n \) components, the
__len__ method is of course useful.
Let us play with some Vec2D objects:
>>> u = Vec2D(0,1)
>>> v = Vec2D(1,0)
>>> w = Vec2D(1,1)
>>> a = u + v
>>> print a
(1, 1)
>>> a == w
True
>>> a = u - v
>>> print a
(-1, 1)
>>> a = u*v
>>> print a
0
>>> print abs(u)
1.0
>>> u == v
False
>>> u != v
True
When you read through this interactive session, you should check
that the calculation is mathematically correct,
that the resulting object type of a calculation is correct,
and how each calculation is performed in the program.
The latter topic is investigated by following the program
flow through the class methods.
As an example, let us consider the expression u != v.
This is a boolean expression that is True since u and v
are different vectors. The resulting object type should be bool,
with values True or False. This is confirmed by the output
in the interactive session above. The Python calculation of
u != v leads to a call to
u.__ne__(v)
which leads to a call to
u.__eq__(v)
The result of this last call is False, because the special method
will evaluate the boolean expression
0 == 1 and 1 == 0
which is obviously False. When going back to the __ne__
method, we end up with a return of not False, which evaluates
to True.
For real computations with vectors in the plane, you would probably just
use a Numerical Python array of length 2. However, one thing such objects
cannot do is evaluating u*v as a scalar product.
The multiplication operator for Numerical Python arrays is not defined
as a scalar product (it is rather defined as
\( (a,b)\cdot(c,d) = (ac, bd) \)). Another difference between our Vec2D
class and Numerical Python arrays is the abs function, which
computes the length of the vector in class Vec2D, while it does something
completely different with Numerical Python arrays.