Basic Usage¶
The most important object in all of kingdon
is Algebra
:
from kingdon import Algebra
An Algebra
has to be initiated with a number of positive, negative,
and null-dimensions, which are traditionally denoted by p
, q
and r
.
For example, in order to create a 2D Geometric Algebra we can initiate
>>> alg = Algebra(p=2, q=0, r=0)
>>> alg = Algebra(2) # Equivalent: default value for p, q, r is 0.
The basis blades of the algebra are available in the dictionary alg.blades
. This can be
added to locals()
in order to allow for easy access to all the basis blades, and allows
the initiation of multivectors using the basis-blades directly:
>>> locals().update(alg.blades)
>>> x = 2 * e + 1 * e1 - 5 * e2 + 6 * e12
where e
is the identity element, i.e. \(e = 1\).
This way of creating multivectors is particularly useful when writing quick scripts
in an interactive environment.
Let’s look at some more general ways of making multivectors, starting with symbolic
multivectors before we go on to numerical multivectors.
Symbolic Multivectors¶
In order to create symbolical multivectors in an algebra, we can call
multivector
and explicitly pass a name
argument.
For example, let us create two symbolic vectors u
and v
in this algebra:
>>> u = alg.multivector(name='u', grades=(1,))
>>> v = alg.multivector(name='v', grades=(1,))
>>> u
(u1) * e1 + (u2) * e2
>>> v
(v1) * e1 + (v2) * e2
The return type of multivector()
is an instance of MultiVector
.
Note
kingdon
offers convenience methods for common types of multivectors, such as the vectors above.
For example, the vectors above can also be created using u = alg.vector(name='u')
.
Moreover, a scalar is created by scalar()
,
a bivector by bivector()
,
a pseudoscalar by pseudoscalar()
, and so on.
However, all of these merely add the corresponding grades
argument to your input and
then call multivector
, so multivector
is what we need to understand.
MultiVector
’s support common math operators:
>>> u + v
(u1 + v1) * e1 + (u2 + v2) * e2
>>> u * v
(u1*v1 + u2*v2) + (u1*v2 - u2*v1) * e12
We also have the inner and exterior “products”:
>>> u | v
(u1*v1 + u2*v2)
>>> u ^ v
(u1*v2 - u2*v1) * e12
We see that in the case of vectors the product is equal to the sum of the inner and exterior.
Since vectors in 2DVGA represent reflections in lines through the origin, we can reflect the
line v
in the line u
by using conjugation:
>>> u >> v
(u1**2*v1 + 2*u1*u2*v2 - u2**2*v1) * e1 + (-u1**2*v2 + 2*u1*u2*v1 + u2**2*v2) * e2
we see that the result is again a vector, as it should be.
These examples should show that the symbolic multivectors of kingdon
make it easy to do symbolic computations. Moreover, we can also use sympy
expressions
as values for the multivector:
>>> from sympy import Symbol, sin, cos
>>> t = Symbol('t')
>>> x = cos(t) * e + sin(t) * e12
>>> x.normsq()
1
More control over basisvectors¶
If we do not just want to create a symbolic multivector of a certain grade,
but with specific blades, we can do so by providing the keys
argument.
>>> x = alg.multivector(name='x', keys=('e1', 'e12'))
>>> (x1) * e1 + (x12) * e12
This can be done either by providing a tuple of strings which indicate which basis-vectors should be present,
or by passing them as integers, i.e. keys=(0b01, 0b11)
is equivalent to the example above.
Internally, kingdon
uses the binary representation.
Numerical Multivectors¶
While kingdon
makes no assumptions about the data structures that are passed into a multivector
in order to support ducktyping and customization as much as possible, it was nonetheless designed to
work really well with numpy
arrays.
For example, to repeat some of the examples above with numerical values, we could do
>>> import numpy as np
>>> uvals, vvals = np.random.random((2, 2))
>>> u = alg.vector(uvals)
>>> v = alg.vector(vvals)
>>> u * v
(0.1541) + (0.0886) * e12
A big performance bottleneck that we suffer from in Python, is that arrays over objects are very slow.
So while we could make a numpy array filled with ~kingdon.multivector.MultiVector
’s, this would tank our performance.
kingdon
gets around this problem by instead accepting numpy arrays as input. So to make a collection of
3 lines, we do
>>> import numpy as np
>>> uvals = np.random.random((2, 3))
>>> u = alg.vector(uvals)
>>> u
([0.82499172 0.71181276 0.98052928]) * e1 + ([0.53395072 0.07312351 0.42464341]) * e2
what is important here is that the first dimension of the array has to have the expected length: 2 for a vector.
All other dimensions are not used by kingdon
. Now we can reflect this multivector in the e1
line:
>>> v = alg.vector((1, 0))
>>> v >> u
([0.82499172 0.71181276 0.98052928]) * e1 + ([-0.53395072 -0.07312351 -0.42464341]) * e2
Despite the different shapes, broadcasting is done correctly in the background thanks to the magic of numpy, and with only minor performance penalties.
Operators¶
Instances of MultiVector
overload all common Geometric Algebra operators.
Below is an overview:
Operation |
Expression |
Infix |
Inline |
---|---|---|---|
Geometric product |
\(ab\) |
|
|
Inner |
\(a \cdot b\) |
|
|
Scalar product |
\(\langle a \cdot b \rangle_0\) |
|
|
Left-contraction |
\(a \rfloor b\) |
|
|
Right-contraction |
\(a \lfloor b\) |
|
|
Outer (Exterior) |
\(a \wedge b\) |
|
|
Regressive |
\(a \vee b\) |
|
|
Conjugate |
\(b a \widetilde{b}\) |
|
|
Project |
\((a \cdot b) \widetilde{b}\) |
|
|
Commutator of |
\(a \times b = \tfrac{1}{2} [a, b]\) |
|
|
Anti-commutator of |
\(\tfrac{1}{2} \{a, b\}\) |
|
|
Sum of |
\(a + b\) |
|
|
Difference of |
\(a - b\) |
|
|
Reverse of |
\(\widetilde{a}\) |
|
|
Squared norm of |
\(a \widetilde{a}\) |
|
|
Norm of |
\(\sqrt{a \widetilde{a}}\) |
|
|
Normalize |
\(a / \sqrt{a \widetilde{a}}\) |
|
|
Square root of |
\(\sqrt{a}\) |
|
|
Dual of |
\(a*\) |
|
|
Undual of |
|
||
Grade |
\(\langle a \rangle_k\) |
|
Note that formally conjugation is defined by \(ba b^{-1}\) and projection by \((a \cdot b) b^{-1}\), but that both are implemented using reversion instead of an inverse. This is because reversion is much faster to calculate, and because in practice \(b\) will often by either a rotor satisfying \(b \widetilde{b} = 1\) or a blade satisfying \(b^2 = b \cdot b\), and thus the inverse is identical to the reverse (up to sign).
If you want to replace these operators by their proper definitions, you can use the register decorator to overwrite the default operator (use at your own risk):
>>> @alg.register(name='sw')
>>> def sw(x, y):
>>> return x * y / y
>>> @alg.register(name='proj')
>>> def proj(x, y):
>>> return (x | y) / y
Graphing using ganja.js
¶
kingdon
supports the ganja.js
graphing syntax. For those already familiar with
ganja.js
, the API will feel very similar:
>>> alg.graph(0xff0000, u, "u", lineWidth=3)
The rules are simple: all positional arguments will be passed on to ganja.js
as
elements to graph, whereas keyword arguments are passed to ganja.js
as options.
Hence, the example above will graph the line u
with lineWidth = 3
,
and will attach the label “u” to it, and all of this will be red.
Identical to ganja.js
, valid inputs to alg.graph
are (lists of) instances
of MultiVector
, strings, and hexadecimal numbers to indicate colors,
or a function without arguments that returns these things.
The strings can be simple labels, or valid SVG syntax.
Note
kingdon supports ganja.js
’s animation and interactivity in jupyter notebooks,
try kingdon in your browser to give it a go!
Performance Tips¶
Because kingdon
attempts to symbolically optimize expressions
using sympy
the first time they are called, the first call to any operation is comparatively slow,
whereas subsequent calls have very good performance.
There are however several things to be aware of to ensure good performance.
Graded¶
The first time kingdon
is asked to perform an operation it hasn’t seen before, it performs code generation
for that particular request. Because codegen is the most expensive step, it is beneficial to reduce the number of
times it is needed. An easy way to achieve this is to initiate the Algebra
with graded=True.
This enforces that kingdon
does not specialize codegen down to the individual basis blades, but rather only
per grade. This means there are far less combinations that have to be considered and generated.
Numba JIT¶
We can enable numba just-in-time compilation by initiating an Algebra
with wrapper=numba.njit.
This comes with a significant cost the first time any operator is called, but subsequent calls to the same operator are
significantly faster. It is worth mentioning that when dealing with Numerical Multivectors over numpy arrays,
the benefit of using numba actually reduces rapidly as the numpy arrays become larger, since then most of the time
is spend in numpy routines anyway.
Register Expressions¶
To make it easy to optimize larger expressions, kingdon
offers the register()
decorator.
>>> alg = Algebra(3, 0, 1)
>>>
>>> @alg.register
>>> def myfunc(u, v):
>>> return u * (u + v)
>>>
>>> x = alg.vector(np.random.random(4))
>>> y = alg.vector(np.random.random(4))
>>> myfunc(x, y)
Calling the decorated myfunc
has the benefit that all the numerical computation is done in one single call,
instead of doing each binary operation individually. This has the benefit that all the (expensive) python boilerplate
code is called only once.