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 𝐞₁ + u2 𝐞₂
>>> v
v1 𝐞₁ + v2 𝐞₂
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) 𝐞₁ + (u2 + v2) 𝐞₂
>>> u * v
(u1*v1 + u2*v2) + (u1*v2 - u2*v1) 𝐞₁₂
We also have the inner and exterior “products”:
>>> u | v
(u1*v1 + u2*v2)
>>> u ^ v
(u1*v2 - u2*v1) 𝐞₁₂
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) 𝐞₁ + (u1**2*v2 - 2*u1*u2*v1 - u2**2*v2) 𝐞₂
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
Note
Strings are also automatically converted to symbolics with SymPy.
So x from the example above can also be created as
>>> x = alg.multivector(e='cos(t)', e12='sin(t)')
>>> 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 𝐞₁ + x12 𝐞₁₂
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) 𝐞₁₂
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]) 𝐞₁ + ([0.53395072 0.07312351 0.42464341]) 𝐞₂
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]) 𝐞₁ + ([-0.53395072 -0.07312351 -0.42464341]) 𝐞₂
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 |
\((-1)^{\text{grade}(b) \text{grade}(a)} 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\) |
|
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!
Large Algebra’s¶
In theory kingdon supports algebra’s up to 36D, but your computer might go up in smoke
if you push it that far. In order to make large’s algebras feasible, kingdon no longer
performs symbolic optimization and caching because this consumes to much memory, and instead
just computes naively.
By default any algebra of \(d > 6\) is considered large, but it can be forced manually with
the large option to Algebra depending on your needs:
>>> alg = Algebra(3, large=True)
>>> alg = Algebra(8, large=False)
For examples of large algebra’s, see the OPNS section of the teahouse.
Performance Tips¶
Because kingdon attempts to symbolically optimize expressions the first time they are called, the first
call to any operation is relatively slow, whereas subsequent calls have extremely good performance.
Note however that even the first execution that includes the symbolic optimization typically
takes on the order of milliseconds, so you probably won’t even notice it.
Since kingdon v2.2.x the symbolic code generation features a port of GAMphetamine.js <https://github.com/enkimute/GAmphetamine.js>_’s very powerful Common Subexpression Elimination (CSE) algorithm, which results in the most optimal code known to man. As in, for those cases in which a hand optimized optimum is known, the CSE optimized code is exactly the same. Moreover, it is quick. Praise Enki.
The table below lists multiplications and additions counted in the emitted Python for representative 3DPGA expressions
when CSE is enabled vs not. Count columns use muls/adds.
Operation |
CSE |
Naive |
|---|---|---|
|
|
|
|
|
|
|
|
|
|
|
|
\(p_1 \vee p_2\) — join of two normalized points (regressive product) |
|
|
\(p \vee \ell\) — normalized point \(p\), line \(\ell\) (bivector) |
|
|
\(p_1 \vee p_2 \vee p_3\) — join of three normalized points |
|
|
\((p \cdot P)\,P^{-1}\) — project normalized point \(p\) on plane \(P\) |
|
|
\((p \cdot P) * P + (p * (1 - P * ~P)).grade(3)\) — normalized plane \(P\) |
|
|
\((p \cdot \ell) * (-\ell) + (p * (1 - \ell * ~ \ell)).grade(3)\) — normalized line \(\ell\) |
|
|
Note
Not all of these counts are achieved yet by kingdon’s default binary operators, but the CSE infrastructure is now in place. In order to achieve these counts out of the box requires to also translate GAmphetamine.js’s type system, which is on the v3.x roadmap.
The symbolically optimized code that kingdon produces is already a good starting point for high performance code. However, there are still several things to be aware of to ensure good performance.
Broadcasting¶
Avoid arrays of multivectors, and use multivectors over e.g. numpy arrays or PyTorch
tensors instead, as shown in Array Syntax.
This ensures the high level overhead of kingdon is paid only once, and we instead delegate
the computation to the underlying datastructures.
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.
Moreover, one can use @alg.register(symbolic=True) to symbolically optimize the expression, similar to how
kingdon’s default binary operators work. As we have seen above in the CSE section, this can result in significant
performance improvements. Afterall, the fastest computation is one you do not have to do.
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 a relatively expensive step, it can be 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,
which will apply numba’s njit decorator to all of kingdon’s generated functions.
This comes with a significant cost the first time any operator is called, but subsequent calls to the same operator are
significantly faster.
However, it is worth mentioning that when dealing with Numerical Multivectors over e.g. numpy arrays,
the benefit of using numba actually disappears rapidly as the numpy arrays become larger, since then most of the time
is spend in numpy routines anyway.
So you need to experiment carefully if numba is right for you.