#!/usr/bin/env python
# -*- coding: utf-8 -*-
"""Mathematical functions and objects.
:Date: 2019-12-12
.. module:: mathematics
:platform: *nix, Windows
:synopsis: Mathematical functions and objects.
.. moduleauthor:: Daniel Weschke <daniel.weschke@directbox.de>
"""
import math
from math import gcd
[docs]def lcm(a, b):
"""Compute the lowest common multiple of a and b"""
return a/gcd(a, b)*b
[docs]class vector(list):
"""Use/create vector like lists
* size, length of a vector use len(a)
* absolute, magntude, norm of a vector use abs(a), see
:meth:`__abs__`
* dot product use a * b, see :meth:`__mul__` and :meth:`__rmul__`
* outer product use a @ b, see :meth:`__matmul__`
:__iter__:
>>> [i*2 for i in vector([1, 2, 3])]
[2, 4, 6]
:__setitem__:
>>> v = vector([1, 2, 3, 4, 5])
>>> v[3:5] = [1, 2]
>>> print(v)
[1, 2, 3, 1, 2]
>>> isinstance(v, vector)
True
:__eq__(a, b):
>>> vector([1, 2, 3, 1, 2]) == vector([1, 2, 3, 1, 2])
True
>>> vector([1, 2, 3, 1, 2]) == vector([1, 2, 3, 2, 1])
False
:__ne__(a, b):
>>> vector([1, 2, 3, 1, 2]) != vector([1, 2, 3, 1, 2])
False
>>> vector([1, 2, 3, 1, 2]) != vector([1, 2, 3, 2, 1])
True
:__contains__:
>>> 2 in vector([1, 2, 3])
True
:__isub__:
>>> v = vector([1, 2, 3])
>>> v -= vector([3, 3, 3])
>>> print(v)
[-2, -1, 0]
:__imul__:
>>> v = vector([1, 2, 3])
>>> v *= vector([3, 3, 3])
>>> print(v)
18
:__imatmul__:
>>> m = vector([1, 2, 3])
>>> m *= vector([3, 3, 3])
>>> print(v)
[[3, 3, 3], [6, 6, 6], [9, 9, 9]]
"""
[docs] def __getitem__(self, index):
"""
For index return value, for range return new vector object.
:Example:
>>> v = vector([1, 2, 3, 4, 5])
>>> v[1:3]
[2, 3]
>>> v = vector([1, 2, 3, 4, 5])
>>> v[3]
4
"""
# use the list.__getitem__ method and convert result to vector
item = super().__getitem__(index)
return vector(item) if isinstance(item, list) else item
[docs] def __pos__(self):
"""\+ a (new object)
.. math::
\mathbf{b} = +\mathbf{a}
"""
return vector([*self])
[docs] def __neg__(self):
"""\- a (new object)
.. math::
\mathbf{b} = -\mathbf{a}
"""
return vector([-i for i in self])
[docs] def __add__(self, other):
"""a + b (new object)
.. math::
\mathbf{c} = \mathbf{a} + \mathbf{b}
:Example:
>>> v = vector([1, 2, 3]) + vector([1, 2, 3])
>>> print(v)
[2, 4, 6]
"""
return vector([i+j for i, j in zip(self, other)])
# overwrite because [1, 2] + [5, 8] = [1, 2, 5, 8]
[docs] def __iadd__(self, other):
"""a += b (new object)
:Example:
>>> v = vector([1, 2, 3])
>>> v += vector([3, 3, 3])
>>> print(v)
[4, 5, 6]
"""
return vector([i+j for i, j in zip(self, other)])
[docs] def __sub__(self, other):
"""a - b (new object)
.. math::
\mathbf{c} = \mathbf{a} - \mathbf{b}
:Example:
>>> v = vector([1, 2, 3]) - vector([1, 2, 3])
>>> print(v)
[0, 0, 0]
"""
return vector([i-j for i, j in zip(self, other)])
[docs] def __mul__(self, other):
r"""Scalar multiplication, dot product (inner product) or
vector-matrix multiplication. (new object)
:type other: scalar, vector (or 1d list), matrix (or 2d list)
.. math::
\mathbf{c} &= \mathbf{a} \cdot b \\
\mathbf{c} &= \mathbf{a} \cdot \mathbf{b} \\
\mathbf{c} &= \mathbf{a} \cdot \mathbf{B}
.. note::
No size checking will be conducted, therefore no exceptions
for wrong usage (result will be nonsense).
:Example:
>>> v = vector([1, 2, 3, 4, 5]) * 3
>>> print(v)
[3, 6, 9, 12, 15]
>>> v = vector([1, 2, 3, 4, 5]) * 3.
>>> print(v)
[3.0, 6.0, 9.0, 12.0, 15.0]
>>> s = vector([1, 2, 3]) * vector([1, 2, 3])
>>> print(s)
14
.. seealso::
:meth:`__rmul__`
"""
try: # vector * vector
return sum([i*j for i, j in zip(self, other)])
except:
try: # vector * matrix
return vector([sum(c*d for c, d in zip(self, b_col)) for
b_col in zip(*other)])
except: # vector * scalar
return vector([i*other for i in self])
[docs] def __rmul__(self, other):
r"""Scalar multiplication, dot product (inner product) or
matrix-vector multiplication. (new object)
:type other: scalar (or 1d list and 2d list)
.. math::
\mathbf{c} &= a \cdot \mathbf{b} \\
\mathbf{c} &= \mathbf{a} \cdot \mathbf{b} \\
\mathbf{c} &= \mathbf{A} \cdot \mathbf{b}
.. note::
No size checking will be conducted, therefore no exceptions
for wrong usage (result will be nonsense).
:Example:
>>> v = 3 * vector([1, 2, 3, 4, 5])
>>> print(v)
[3, 6, 9, 12, 15]
>>> v = 3. * vector([1, 2, 3, 4, 5])
>>> print(v)
[3.0, 6.0, 9.0, 12.0, 15.0]
.. seealso::
:meth:`__mul__` and :meth:`matrix.__mul__` for matrix * vector
"""
try: # 2d list * vector (matrix * vector see matrix.__mul__)
return vector([sum(c*d for c, d in zip(a_row, self)) for
a_row in other])
except: # scalar * vector
return self*other
[docs] def __matmul__(self, other):
r"""Outer product a @ b (new object)
.. math::
\mathbf{c} = \mathbf{a} \otimes \mathbf{b}
=
\begin{pmatrix}
a_{1}\\
a_{2}\\
a_{3}
\end{pmatrix}
\otimes
\begin{pmatrix}
b_{1}\\
b_{2}\\
b_{3}
\end{pmatrix}
=
\begin{pmatrix}
a_{1}b_{1}&a_{1}b_{2}&a_{1}b_{3}\\
a_{2}b_{1}&a_{2}b_{2}&a_{2}b_{3}\\
a_{3}b_{1}&a_{3}b_{2}&a_{3}b_{3}
\end{pmatrix}
:Example:
>>> m = vector([1, 2, 3]) @ vector([1, 2, 4])
>>> print(m)
[[1, 2, 4], [2, 4, 8], [3, 6, 12]]
"""
try: # vector * vector
return matrix([[i*j for j in other] for i in self])
except: # vector * number
return self*other
[docs] def __abs__(self):
r"""Magnitude / norm of a vector
.. math::
b = |\mathbf{a}| = \sqrt{\sum a_i^2} =
\sqrt{\mathbf{a} \cdot \mathbf{a}}
:Example:
>>> v = vector([3, 4])
>>> abs(v)
5.0
"""
return math.sqrt(self * self)
[docs] def __lt__(self, other):
"""Test if this object is lower (smaller) than the other object.
.. math::
|\mathbf{a}| \lt |\mathbf{b}|
:Example:
>>> vector([3, 2, 1]) < vector([1, 2, 3])
False
>>> vector([3, 2, 1]) < vector([1, 2, 4])
True
"""
return abs(self) < abs(other)
[docs] def __le__(self, other):
"""Test if this object is lower (smaller) than or equal the
other object.
.. math::
|\mathbf{a}| \le |\mathbf{b}|
:Example:
>>> vector([3, 2, 1]) <= vector([1, 2, 3])
True
>>> vector([3, 2, 1]) <= vector([1, 2, 2])
False
"""
return abs(self) <= abs(other)
[docs] def __gt__(self, other):
"""Test if this object is greater (larger) than the other
object.
.. math::
|\mathbf{a}| \gt |\mathbf{b}|
:Example:
>>> vector([1, 2, 3]) > vector([3, 2, 1])
False
>>> vector([1, 2, 3]) > vector([2, 2, 1])
True
"""
return abs(self) > abs(other)
[docs] def __ge__(self, other):
"""Test if this object is greater (larger) than or equal the
other object.
.. math::
|\mathbf{a}| \ge |\mathbf{b}|
:Example:
>>> vector([1, 2, 3]) >= vector([3, 2, 1])
True
>>> vector([1, 2, 3]) >= vector([4, 2, 1])
False
"""
return abs(self) >= abs(other)
[docs] def __str__(self):
return str([*self])
[docs] def __repr__(self):
return "vector(" + str(self) + ")"
[docs] @staticmethod
def full(length, fill_value):
"""Returns a vector of length m filled with v.
:param length: length of the vector, e. g. 3
:type length: int
:param fill_value: fill value
:Type fill_value: scalar
:Example:
>>> v = vector.full(3, 7)
>>> print(v)
[7, 7, 7]
"""
return vector([fill_value for i in range(length)])
[docs] @staticmethod
def zeros(length):
"""Returns a zero vector of length m filled with zeros.
:param length: length of the vector, e. g. 3
:type length: int
:Example:
>>> v = vector.zeros(3)
>>> print(v)
[0.0, 0.0, 0.0]
"""
return vector.full(length, 0.)
[docs] @staticmethod
def ones(length):
"""Returns a vector of length m filled with ones.
:param length: lhape of the vector, e. g. 3
:type length: int
:Example:
>>> v = vector.ones(3)
>>> print(v)
[1.0, 1.0, 1.0]
"""
return vector.full(length, 1.)
[docs] @staticmethod
def random(length, lmin=0.0, lmax=1.0):
"""Returns a random vector of length n filled with random
numbers.
:param length: lhape of the vector, e. g. 3
:type length: int
:Example:
>>> v = vector.random(3)
>>> print(v)
[0.9172905912930438, 0.8908124278322492, 0.5256002790725927]
>>> v = vector.random(3, 1, 2)
>>> print(v)
[1.2563665665080803, 1.9270454509964547, 1.2381672401270487]
"""
import random
dl = lmax-lmin
return vector([dl*random.random()+lmin for i in range(length)])
[docs] @staticmethod
def normalized(a):
r"""Normalize a vector (i. e. the vector has a length of 1)
:type a: vector
.. math::
\mathbf{\hat{a}} = \frac{\mathbf{a}}{|\mathbf{a}|}
.. seealso::
:meth:`__abs__` for a norm (magnitude) of a vector
"""
a_mag = abs(a)
return vector([i/a_mag for i in a])
[docs] @staticmethod
def ang(a, b):
x = a * b / (abs(a) * abs(b))
# decimal floating-point numbers are only approximated by the
# binary floating-point numbers actually stored in the machine.
# https://docs.python.org/3.8/tutorial/floatingpoint.html
xr15 = round(x, 15)
return math.acos(xr15) if -1 <= xr15 <= 1 else math.acos(x)
[docs] @staticmethod
def conjugate(a):
"""
New vector object
:type a: list
"""
return vector([i.conjugate() for i in a])
[docs] @staticmethod
def re(a):
"""Return the real parts of a complex vector
:type a: list
"""
return vector([i.real for i in a])
[docs] @staticmethod
def im(a):
"""Return the imaginary parts of a complex vector
:type a: list
"""
return vector([i.imag for i in a])
[docs] @staticmethod
def abs(a):
"""Return modulus parts of a complex vector
:type a: list
"""
return vector([abs(i) for i in a])
[docs] @staticmethod
def arg(a):
"""Return phase parts of a complex vector
:type a: list
"""
return vector([math.atan2(i.imag, i.real) for i in a])
[docs] @staticmethod
def cross(a, b):
r"""Cross product
:type a: list
:type b: list
c is orthogonal to both a and b.
The direction of c can be found with the right-hand rule.
.. math::
\mathbf{c} = \mathbf{a} \times \mathbf{b}
"""
return vector([a[1]*b[2] - a[2]*b[1],
a[2]*b[0] - a[0]*b[2],
a[0]*b[1] - a[1]*b[0]])
[docs] @staticmethod
def isclose(a, b, rel_tol=0.05, abs_tol=1e-8):
return all([math.isclose(i, j, rel_tol=rel_tol, abs_tol=abs_tol) for i, j in zip(a, b)])
[docs] def iscloseto(self, other):
return vector.isclose(self, other)
[docs] def normalize(self):
r"""Normalize a vector (i. e. the vector has a length of 1)
:type a: vector
.. math::
\mathbf{\hat{a}} = \frac{\mathbf{a}}{|\mathbf{a}|}
.. seealso::
:meth:`__abs__` for a norm (magnitude) of a vector
"""
self[:] = vector.normalized(self)
return self
[docs] def rotate_x(self, theta):
r"""Rotation about the x dirction.
.. math::
\begin{bmatrix}x' \\ y' \\ z' \\ h'\end{bmatrix} =
\begin{bmatrix}
1 & 0 & 0 & 0 \\
0 & \cos \theta & -\sin \theta & 0 \\
0 & \sin \theta & \cos \theta & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}x \\ y \\ z \\ h\end{bmatrix}
"""
self[:] = matrix.rx(theta) * self
return self
[docs] def rotate_y(self, theta):
r"""Rotation about the y dirction.
.. math::
\begin{bmatrix}x' \\ y' \\ z' \\ h'\end{bmatrix} =
\begin{bmatrix}
\cos \theta & 0 & \sin \theta & 0 \\
0 & 1 & 0 & 0 \\
-\sin \theta & 0 & \cos \theta & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}x \\ y \\ z \\ h\end{bmatrix}
"""
self[:] = matrix.ry(theta) * self
return self
[docs] def rotate_z(self, theta):
r"""Rotation about the z dirction.
.. math::
\begin{bmatrix}x' \\ y' \\ z' \\ h'\end{bmatrix} =
\begin{bmatrix}
\cos \theta & -\sin \theta & 0 & 0 \\
\sin \theta & \cos \theta & 0 & 0 \\
0 & 0 & 1 & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}x \\ y \\ z \\ h\end{bmatrix}
"""
self[:] = matrix.rz(theta) * self
return self
[docs] def translate(self, tx, ty, tz):
r"""Translation
.. math::
\begin{bmatrix}x' \\ y' \\ z' \\ h'\end{bmatrix} =
\begin{bmatrix}
1 & 0 & 0 & t_x \\
0 & 1 & 0 & t_y \\
0 & 0 & 1 & t_z \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}x \\ y \\ z \\ h\end{bmatrix}
"""
self[:] = matrix.t(tx, ty, tz) * self
return self
[docs] def scale(self, sx, sy=None, sz=None):
r"""Scaling
uniform scaling if sx=sy=sz=s.
Note that scaling happens around the origin, so objects not
centered at the origin will have their centers move. To avoid
this, either scale the object when it's located at the origin,
or perform a translation afterwards to move the object back to
where it should be.
.. math::
\begin{bmatrix}x' \\ y' \\ z' \\ h'\end{bmatrix} =
\begin{bmatrix}
s_x & 0 & 0 & 0 \\
0 & s_y & 0 & 0 \\
0 & 0 & s_z & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}x \\ y \\ z \\ h\end{bmatrix}
"""
# if not sy is not suitable because 0 is also false
if sy is None:
sy = sx
sz = sx
self[:] = matrix.s(sx, sy, sz) * self
return self
[docs] def ch_cs(self, cs):
r"""Transform this vector from its defined coordinate system to
a new coordinate system, defined by the given coordinate system
(u, v and w direction vectors).
.. math::
\begin{bmatrix}x' \\ y' \\ z' \\ h'\end{bmatrix} =
\begin{bmatrix}
u_x & u_y & u_z & 0 \\
v_x & v_y & v_z & 0 \\
w_x & w_y & w_z & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}x \\ y \\ z \\ h\end{bmatrix}
"""
self[:] = cs * self
return self
[docs]class matrix(list):
"""Use/create matrix like list of lists
"""
[docs] def __getitem__(self, index):
r"""
For index return value, for range return new vector object.
:Example:
>>> m = matrix([[1, 2, 3, 0], [4, 5, 6, 0], [7, 8, 9, 0], \
[0, 0, 0, 0]])
>>> print(m[:])
[[1, 2, 3, 0], [4, 5, 6, 0], [7, 8, 9, 0], [0, 0, 0, 0]]
>>> print(m[2])
[7, 8, 9, 0]
>>> print(m[2, :])
[7, 8, 9, 0]
>>> print(m[:, 2])
[3, 6, 9, 0]
>>> print(m[2, 2])
9
>>> print(m[:, 1:3])
[[2, 3], [5, 6], [8, 9], [0, 0]]
>>> print(m[0:2, 1:3])
[[2, 3], [5, 6]]
>>> print(m[::2, ::2])
[[1, 3], [7, 9]]
"""
# TODO single row or column = vector (1d list)?
# index: slice(stop), slice(start, stop[, step])
# for m[(1,3),:] -> index = ((1, 3), slice(None, None, None))
# use the list.__getitem__ method and convert result to matrix
#print(index)
try: # 2d slicing (tuple of sclices)
try: # range:range or range:single: [:, 2]
#print(1)
item = [row.__getitem__(index[1]) for
row in super().__getitem__(index[0])]
except: # single:range: [2, :], [2, 2]
#print(2)
item = super().__getitem__(index[0]).__getitem__(index[1])
except: # 1d slicing: [:], [2], [2][2]
#print(3)
item = super().__getitem__(index)
#print(item)
return matrix(item) if isinstance(item, list) and \
isinstance(item[0], list) else item
[docs] def __setitem__(self, index, value):
#print(index, value)
try: # [2, 2]
super().__getitem__(index[0]).__setitem__(index[1], value)
except: # [2], [2][2]
super().__setitem__(index, value)
#print(self)
[docs] def __mul__(self, other):
r"""Scalar multiplication, dot product (inner product) or
matrix-vector multiplication. (new object)
:type other: scalar, vector (or 1d list), matrix (or 2d list)
.. math::
\mathbf{C} &= \mathbf{A} \cdot b \\
\mathbf{c} &= \mathbf{A} \cdot \mathbf{b} \\
\mathbf{C} &= \mathbf{A} \cdot \mathbf{B}
.. note::
No size checking will be conducted, therefore no exceptions
for wrong usage (result will be nonsense).
:Example:
>>> m = matrix([[1, 2, 3, 4], [5, 6, 7, 8], [9, 10, 11, 12], \
[0, 0, 0, 1]]) * 5
>>> print(m)
[[5, 10, 15, 20], [25, 30, 35, 40], [45, 50, 55, 60], [0, 0, 0, 5]]
>>> m = matrix([[1, 2, 3, 4], [5, 6, 7, 8], [9, 10, 11, 12], \
[0, 0, 0, 1]]) * 5.
>>> print(m)
[[5.0, 10.0, 15.0, 20.0], [25.0, 30.0, 35.0, 40.0], [45.0, 50.0, 55.0, 60.0], [0.0, 0.0, 0.0, 5.0]]
>>> v = matrix([[1, 2, 3], [4, 5, 6], [7, 8, 9]]) * \
vector([12, 12, 13])
>>> print(v)
[75, 186, 297]
>>> m = matrix([[1, 2, 3], [4, 5, 6], [7, 8, 9]]) * \
matrix([[12, 12, 13], [14, 15, 16], [17, 18, 19]])
>>> print(m)
[[91, 96, 102], [220, 231, 246], [349, 366, 390]]
.. seealso::
:meth:`__rmul__`
"""
try: # matrix * matrix
return matrix([[sum(c*d for c, d in zip(a_row, b_col)) for
b_col in zip(*other)] for a_row in self])
except:
try: # matrix * vector
return vector([sum(c*d for c, d in zip(a_row, other)) for
a_row in self])
except: # matrix * scalar
return matrix([[a*other for a in a_row] for a_row in self])
[docs] def __rmul__(self, other):
r"""Scalar multiplication, dot product (inner product) or
vector-matrix multiplication. (new object)
:type other: scalar (or 1d list and 2d list)
.. math::
\mathbf{C} &= a \cdot \mathbf{B} \\
\mathbf{c} &= \mathbf{a} \cdot \mathbf{B} \\
\mathbf{C} &= \mathbf{A} \cdot \mathbf{B}
.. note::
No size checking will be conducted, therefore no exceptions
for wrong usage (result will be nonsense).
:Example:
>>> m = 5 * matrix([[1, 2, 3, 4], [5, 6, 7, 8], \
[9, 10, 11, 12], [0, 0, 0, 1]])
>>> print(m)
[[5, 10, 15, 20], [25, 30, 35, 40], [45, 50, 55, 60], [0, 0, 0, 5]]
>>> m = 5. * matrix([[1, 2, 3, 4], [5, 6, 7, 8], \
[9, 10, 11, 12], [0, 0, 0, 1]])
>>> print(m)
[[5.0, 10.0, 15.0, 20.0], [25.0, 30.0, 35.0, 40.0], [45.0, 50.0, 55.0, 60.0], [0.0, 0.0, 0.0, 5.0]]
.. seealso::
:meth:`__mul__` and :meth:`vector.__mul__` for vector * matrix
"""
try: # 2d list * matrix (matrix * matrix see matrix.__mul__)
return matrix([[sum(c*d for c, d in zip(a_row, b_col)) for
b_col in zip(*self)] for a_row in other])
except:
try: # 1d list * matrix (vector * matrix see vector.__mul__)
return vector([sum(c*d for c, d in zip(other, b_col)) for
b_col in zip(*self)])
except: # scalar * vector
return self*other
[docs] def __str__(self):
return str([*self])
[docs] def __repr__(self):
return "matrix(" + str(self) + ")"
[docs] @staticmethod
def zeros(m, n):
"""Returns a zero matrix of size mxn; m rows and n columns
filled with zeros.
:param m: number of rows of the matrix, e. g. 3
:type m: int
:param n: number of columns of the matrix, e. g. 3
:type n: int
:Example:
>>> m = matrix.zeros(3, 3)
>>> print(m)
[[0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]]
"""
return matrix([[0. for i in range(n)] for j in range(m)])
[docs] @staticmethod
def rx(theta):
r"""Rotation matrix about the x direction.
Rotates the coordinate system of vectors
.. math::
\mathbf{R}_{x}(\theta) =
\begin{bmatrix}
1 & 0 & 0 & 0 \\
0 & \cos \theta & -\sin \theta & 0 \\
0 & \sin \theta & \cos \theta & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
"""
s = math.sin(theta)
c = math.cos(theta)
T = matrix([[ 1, 0, 0, 0],
[ 0, c,-s, 0],
[ 0, s, c, 0],
[ 0, 0, 0, 1]])
return T
[docs] @staticmethod
def ry(theta):
r"""Rotation matrix about the y direction.
Rotates the coordinate system of vectors
.. math::
\mathbf{R}_{y}(\theta) =
\begin{bmatrix}
\cos \theta & 0 & \sin \theta & 0 \\
0 & 1 & 0 & 0 \\
-\sin \theta & 0 & \cos \theta & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
"""
s = math.sin(theta)
c = math.cos(theta)
T = matrix([[ c, 0, s, 0],
[ 0, 1, 0, 0],
[-s, 0, c, 0],
[ 0, 0, 0, 1]])
return T
[docs] @staticmethod
def rz(theta):
r"""Rotation matrix about the z direction.
Rotates the coordinate system of vectors
.. math::
\mathbf{R}_{z}(\theta) =
\begin{bmatrix}
\cos \theta & -\sin \theta & 0 & 0 \\
\sin \theta & \cos \theta & 0 & 0 \\
0 & 0 & 1 & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
:Example:
.. math::
\mathbf{R}_{z}(\theta)
\begin{bmatrix}1 \\ 0 \\ 0 \\ 1\end{bmatrix}
&=
\begin{bmatrix}
\cos 90° & -\sin 90° & 0 & 0 \\
\sin 90° & \cos 90° & 0 & 0 \\
0 & 0 & 1 & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}1 \\ 0 \\ 0 \\ 1\end{bmatrix} \\
&=
\begin{bmatrix}
0 & -1 & 0 & 0 \\
1 & 0 & 0 & 0 \\
0 & 0 & 1 & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}1 \\ 0 \\ 0 \\ 1\end{bmatrix} =
\begin{bmatrix}0 \\ 1 \\ 0 \\ 1\end{bmatrix}
"""
s = math.sin(theta)
c = math.cos(theta)
T = matrix([[ c,-s, 0, 0],
[ s, c, 0, 0],
[ 0, 0, 1, 0],
[ 0, 0, 0, 1]])
return T
[docs] @staticmethod
def t(tx, ty, tz):
r"""Translation matrix
.. math::
\mathbf{T} =
\begin{bmatrix}
1 & 0 & 0 & t_x \\
0 & 1 & 0 & t_y \\
0 & 0 & 1 & t_z \\
0 & 0 & 0 & 1
\end{bmatrix}
"""
T = matrix([[ 1, 0, 0,tx],
[ 0, 1, 0,ty],
[ 0, 0, 1,tz],
[ 0, 0, 0, 1]])
return T
[docs] @staticmethod
def s(sx, sy=None, sz=None):
r"""Scaling matrix
uniform scaling if sx=sy=sz=s.
Note that scaling happens around the origin, so objects not
centered at the origin will have their centers move. To avoid
this, either scale the object when it's located at the origin,
or perform a translation afterwards to move the object back to
where it should be.
.. math::
\mathbf{S} =
\begin{bmatrix}
s_x & 0 & 0 & 0 \\
0 & s_y & 0 & 0 \\
0 & 0 & s_z & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
"""
# if not sy is not suitable because 0 is also false
if sy is None:
sy = sx
sz = sx
T = matrix([[sx, 0, 0, 0],
[ 0,sy, 0, 0],
[ 0, 0,sz, 0],
[ 0, 0, 0, 1]])
return T
[docs] @staticmethod
def transposed(a):
"""Transpose
"""
m = len(a)
n = len(a[0])
new = matrix.zeros(n, m)
for i in range(m):
for j in range(n):
new[j][i] = a[i][j]
return new
[docs] def transpose(self):
"""Transpose
"""
self[:] = matrix.transposed(self)
return self
[docs] def rotate_x(self, theta):
r"""Rotation about the x dirction.
.. math::
\begin{bmatrix}
xx & xy & xz & t_x \\
yx' & yy' & yz' & t_y' \\
zx' & zy' & zz' & t_z' \\
0 & 0 & 0 & h
\end{bmatrix} =
\begin{bmatrix}
1 & 0 & 0 & 0 \\
0 & \cos \theta & -\sin \theta & 0 \\
0 & \sin \theta & \cos \theta & 0 \\
0 & 0 & 0 & h
\end{bmatrix}
\begin{bmatrix}
xx & xy & xz & t_x \\
yx & yy & yz & t_y \\
zx & zy & zz & t_z \\
0 & 0 & 0 & h
\end{bmatrix}
"""
self[:] = matrix.rx(theta) * self
return self
[docs] def rotate_y(self, theta):
r"""Rotation about the y dirction.
.. math::
\begin{bmatrix}
xx' & xy' & xz' & t_x' \\
yx & yy & yz & t_y \\
zx' & zy' & zz' & t_z' \\
0 & 0 & 0 & h
\end{bmatrix} =
\begin{bmatrix}
\cos \theta & 0 & \sin \theta & 0 \\
0 & 1 & 0 & 0 \\
-\sin \theta & 0 & \cos \theta & 0 \\
0 & 0 & 0 & h
\end{bmatrix}
\begin{bmatrix}
xx & xy & xz & t_x \\
yx & yy & yz & t_y \\
zx & zy & zz & t_z \\
0 & 0 & 0 & h
\end{bmatrix}
"""
self[:] = matrix.ry(theta) * self
return self
[docs] def rotate_z(self, theta):
r"""Rotation about the z dirction.
.. math::
\begin{bmatrix}
xx' & xy' & xz' & t_x' \\
yx' & yy' & yz' & t_y' \\
zx & zy & zz & t_z \\
0 & 0 & 0 & h
\end{bmatrix} =
\begin{bmatrix}
\cos \theta & -\sin \theta & 0 & 0 \\
\sin \theta & \cos \theta & 0 & 0 \\
0 & 0 & 1 & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}
xx & xy & xz & t_x \\
yx & yy & yz & t_y \\
zx & zy & zz & t_z \\
0 & 0 & 0 & h
\end{bmatrix}
"""
self[:] = matrix.rz(theta) * self
return self
[docs] def translate(self, tx, ty, tz):
r"""Translation
.. math::
\begin{bmatrix}
xx & xy & xz & t_x' \\
yx & yy & yz & t_y' \\
zx & zy & zz & t_z' \\
0 & 0 & 0 & h
\end{bmatrix} =
\begin{bmatrix}
1 & 0 & 0 & t_x \\
0 & 1 & 0 & t_y \\
0 & 0 & 1 & t_z \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}
xx & xy & xz & t_x \\
yx & yy & yz & t_y \\
zx & zy & zz & t_z \\
0 & 0 & 0 & h
\end{bmatrix}
"""
self[:] = matrix.t(tx, ty, tz) * self
return self
[docs] def scale(self, sx, sy=None, sz=None):
r"""Scaling
uniform scaling if sx=sy=sz=s.
Note that scaling happens around the origin, so objects not
centered at the origin will have their centers move. To avoid
this, either scale the object when it's located at the origin,
or perform a translation afterwards to move the object back to
where it should be.
.. math::
\begin{bmatrix}
xx' & xy' & xz' & t_x' \\
yx' & yy' & yz' & t_y' \\
zx' & zy' & zz' & t_z' \\
0 & 0 & 0 & h
\end{bmatrix} =
\begin{bmatrix}
s_x & 0 & 0 & 0 \\
0 & s_y & 0 & 0 \\
0 & 0 & s_z & 0 \\
0 & 0 & 0 & 1
\end{bmatrix}
\begin{bmatrix}
xx & xy & xz & t_x \\
yx & yy & yz & t_y \\
zx & zy & zz & t_z \\
0 & 0 & 0 & h
\end{bmatrix}
"""
# if not sy is not suitable because 0 is also false
if sy is None:
sy = sx
sz = sx
self[:] = matrix.s(sx, sy, sz) * self
return self