Geometric Algebra Library

10 Jun 2004

      Hi,

[crossposted to oo-numerics and boost-developers]

About 20 months back, I asked if there was any interest in a geometric or
Clifford algebra library written in C++. Despite some positive replies, I
got the
impression that there were not a lot of people in the Boost community who
saw any practical uses for such library.

Back then I had a rudimentary proof-of-concept for a Clifford algebra
library that used meta-programming techniques to gain significant
performance benefits over other existing libraries (Gaigen, nKlein, CLU,
etcetera).

Over the past few months I have put a lot more work into this library. The
bulk of the work consisted of adding more unit-tests, proving many of the
theorems in some of the papers by David Hestenes et al. I also added
expression templates to improve performance.

To give an idea of what the library is capable of:

#include <boost/math/clifford/clifford.hpp>
#include <boost/mpl/list_c.hpp>
#include <cassert>

using namespace boost::math::clifford;
using namespace boost::mpl;

typedef algebra<4, 1> conformal_model;

typedef blade<float, 0, conformal_model> scalar;
typedef blade<float, 1, conformal_model> vector;
typedef blade<float, 5, conformal_model> pseudo_scalar;
typedef blade<float, 2, conformal_model> bi_vector;

// same as:
// typedef element<float, list_c<int, 1>, conformal_model> bi_vector;

typedef list_c<int, 0, 2> spinor_grades;
typedef element<float, spinor_grades, conformal_model> spinor;

typedef list_c<int, 0, 1, 2, 3, 4, 5> multi_vector_grades;
typedef element<float, multi_vector_grades, conformal_model> multi_vector;

int main()
{
    assert(sizeof(scalar) == sizeof(float));
    assert(sizeof(vector) == 5 * sizeof(float));
    assert(sizeof(bi_vector) == 10 * sizeof(float));
    assert(sizeof(pseudo_scalar) == sizeof(float));
    assert(sizeof(spinor) == 11 * sizeof(float));
    assert(sizeof(multi_vector) == 32 * sizeof(float));

    vector e_plus;
    e_plus[3] = 1.0;

    vector e_min;
    e_min[4] = 1.0;

    e_star = e_min - e_plus;
    e_star *= -std::sqrt(2.0);
    e_star /= 2.0f;
    e [noalias] = e_plus * std::sqrt(2.0) - e_star;

    bi_vector E = outer_product(e, e_star);

    // geometric product between bi_vector and vector
    multi_vector v = E * e_star;

    // geometric product between two vectors, generating a spinor
    vector a, b;
    spinor s = a * b;

    // doing a spinor rotation
    vector c;
    c [noalias] = reverse(s) * a * s;
}

Each product always does the bare minimum of additions, multiplications and
subtractions to calculate the result on the left-hand-side of an expression.
Thus, if you do the full geometric product of two multi_vectors but assign
it to a scalar, it will only do those calculations that are involved in
grade 0.

Also, unlike other existing libraries, no dynamic allocations are being done
ever.

Currently, the following operations are supported, between any two random
elements of an algebra (blades and composite multi_vectors):

    geometric_product = operator *
    inner_product
    outer_product
    fat_dot_product
    left_contraction
    right_contraction
    scalar_product

of course, all the normal operations (scalar multiplication, addition of
elements, etcetera) are supported, as well as reversion, involution and
conjugates.

Adding more operators is trivial. They can either be implemented in terms of
existing operators, or by implementing an appropriate grade-predicate. For
example, this is what the implementation of the outer_product looks like:

    struct outer_product_predicate
    {
        template <class ResultGrade, class LhsGrade, class RhsGrade>
        struct apply
        {
            typedef typename mpl::if_
            <
                typename mpl::equal_to
                <
                    ResultGrade,
                    typename mpl::plus<LhsGrade,RhsGrade>::type
                >::type,
                mpl::true_,
                mpl::false_
            >::type type;
        };
    };

    BOOST_MATH_CLIFFORD_PRODUCT_IMPL(
        outer_product,
        impl::outer_product_predicate
    )

which for somebody that knows the definition of the outer_product is easy to
understand.

I have unit-tests for 2, 3 4 and 5 dimensional algebras, but more is
possible if your compiler supports it.

Currently I have only tested the library on MSVC 7.1. I tried 7.0 but
couldn't get it to work (excessive use of MPL). GCC and Intel compiled it at
one point, but because I only use MSVC I dropped support. Of course, if
anybody is interested I'd be trying other compilers again. You will need a
compiler that does aggresive inlining though. On MSVC I am using the
force_inline feature to get the release performance I need.

Anyway, let me know if anybody is interested, and I'll put together a
package to upload somewhere.

Cheers,

Jaap

Jaap Suter

Ken Hagan

Hubert Holin

Petr Kočmíd

Justinas V.D.

Paul C. Leopardi

David B. Held

tags

participants (7)