GluCat: Generic library of universal Clifford algebra templates
| |
||
GluCat is a library of template classes which model the universal Clifford algebras over the field of real numbers, with arbitrary dimension and arbitrary signature. GluCat implements a model of each Clifford algebra corresponding to each non-degenerate quadratic form up to a maximum number of dimensions.
GluCat classes were originally designed to be used as template parameters for other template libraries, such as Blitz++, deal.II, Matrix Template Library and POOMA. These template libraries expect a numeric class which implements "appropriate numeric semantics". To provide these semantics, the GluCat interface matches float or complex<> as much as possible, but the GluCat classes must be different from float or complex<>, because they model a different algebra.
The main template classes are index_set<LO,HI>, framed_multi<Scalar_T,LO,HI> and matrix_multi<Scalar_T,LO,HI>. Scalar_T is a class which represents a number field, such as float or double. LO, an Index_T <= 0 and HI, an Index_T >= 0, are the lower and upper limits for indices in index sets. Testing so far has been done with Index_T as int and Scalar_T as float, double and long double.
The notation for index_set<LO,HI> looks like: {-2,3,4}, a set of integers, with templated limits. 0 is never a member of the set. Each index set represents the non-commutative product of the generators corresponding to each index in strictly increasing order. This is usually called a strictly increasing multi-index in the literature.
eg. {-2,3,4} == {-2}*{3}*{4}
Input treats index sets as sets, so {-2,3,4} == {4,3,-2}. Output is in increasing order.
A framedmultivector, a multivector expressed with respect to a basis, "knows" which subalgebra contains it. It is the subalgebra defined by the frame. A frame here is an ordered set of basis vectors which generates a subalgebra.
The notation for framed_multi<> and matrix_multi<> is the framed notation: 1+2{-3}+4{-2}+5{-3,-2,2}, ie a sum of terms, where each term consists of a coordinate and an index_set.
framed_multi<float,LO,HI> and matrix_multi<float,LO,HI> represent elements of subalgebras of R(HI,-LO).
{-2}*{-2} == -1. {2}*{2} == 1.
The framed notation automatically takes care of the signature of the generators. You don't need to specify a signature to use framed_multi<> or matrix_multi<>, ie. the value of a multivector in this notation is unambiguous.
The framed notation, with negative indices, is not new. Mark Ashdown uses negative indices in GA Package for Maple V. G. P. Wene uses negative indices in his paper on infinite dimensional Clifford algebras. Leo Dorst defines various products in terms of index sets in his paper on "Honing geometric algebras". Finally, and most importantly, Arvind Raja's design and prototype code (the basis for GluCat) uses index sets internally.
framed_multi<> and matrix_multi<> provide the usual arithmetic operators +, -, * and /, as well as % for left contraction, & for Hestenes inner product and ^ for outer product. They also provide operator[] for subscripting using an index set, and operator() for grading.
framed_multi<> and matrix_multi<> also provide the usual functions, conj, real, imag, abs, norm, pow, exp, sqrt, log, cos, acos, cosh, acosh, sin, asin, sinh, asinh, tan, atan, tanh, atanh, plus multivector functions even, frame, inv, involute, reverse, pure, quad, scalar etc.
The aim is to get the notation as close as possible to that used by float or complex<>, while not straying too far from current notations for Clifford algebras.
One of the beauties of the Clifford algebras is that F(P,Q) contains subalgebras isomorphic to F(p,q) for p<=P and q<=Q.
It is therefore OK to add and multiply multivectors from different subalgebras. The result is, in general, contained in the subalgebra which spans the operands.
eg. (2{-1}+3{-2}) + (4{1,2}+5{1,4}) == 2{-1}+3{-2}+4{1,2}+5{1,4}
The idea in using the GluCat classes is to make the range LO...HI as large as you will need for your program. You can use the same class, say framed_multi<double,-4,4> to perform calculations in R(0,2) (quaternions), R(3,0) (3D geometric algebra), R(0,3) (3D Clifford algebra used in Clifford analysis) and R(1,3) (spacetime algebra).
Even better, you can just let LO and HI take default values, which are large enough for most purposes. The defaults are based on sizeof(set_value_t), where set_value_t is defined as unsigned long. On the Intel i386 architecture, unsigned long is 32 bits, DEFAULT_LO is -16 and DEFAULT_HI is 16. So to accommodate the subalgebras listed above, you can just use framed_multi<double>. This incurs minor overhead in size and speed when compared to framed_multi<double,-4,4>. On AMD 64 or Compaq Alpha, unsigned long is 64 bits, DEFAULT_LO is -32 and DEFAULT_HI is 32.
Conversion between framed_multi<> and matrix_multi<> is necessary in many circumstances. matrix_multi<> cannot operate directly on operands in different subalgebras, so such operations work by converting the operands to a common subalgebra by conversion to and from framed_multi<>.
Past a threshold number of dimensions, this conversion is implemented using a generalized Fast Fourier Transform.
GluCat adapts Arvind Raja's prototype code, with permission, for the basic operations and adds /, inv, sqrt and the transcendental functions listed above.
The class index_set<LO,HI> is based on std::bitset<HI-LO+1>.
The class framed_multi<Scalar_T,LO,HI> is based on std::map<index_set<LO,HI>,Scalar_T> or (optionally) on __gnu_cxx::hash_map< const index_set<LO,HI>,Scalar_T, hash <LO,HI> > or std::tr1::unordered_map< const index_set<LO,HI>,Scalar_T, hash <LO,HI> > >. An object of class framed_multi<> represents a multivector as a sum of terms, where each term consists of a basis element and a coordinate. You can also think of this as a map from the set of index sets to the set of coordinate values. Where necessary, framed_multi<> implementation uses std::pair<index_set<LO,HI>,Scalar_T> to represent terms.
The class matrix_multi<> consists of a matrix, plus an index_set<> which represents the frame which defines the subalgebra. matrix_multi<> is currently based on the ublas:: templates provided by the Boost uBLAS Library. GluCat currently uses either sparse matrices in row compressed format or dense matrices.
The class matrix_multi<> uses the matrix generators of the minimum real representation of R(p,q) as per Porteous and Lounesto. A single instance of the class generator_table<> is used to construct and store the generators. This class, based on std::map<>, maps from signatures (p,q) to vectors of matrix generators.
Division and inverse currently use LU decomposition, plus iterative refinement, as per Higham. The class matrix_multi<> uses framed_multi<> for most operations other than +, -, *, /.
Conversion from framed_multi<> to matrix_multi<> uses the smallest subalgebra which contains the framed_multi<>.
The fold() and unfold() operations map an index set within a subalgebra to and from the isomorphic subalgebra defined by the corresponding contiguous index set. An index set defines a subalgebra by defining a frame, which is an ordered set of basis vectors. The indices in the index set, in order, define this ordered set. In turn, the ordered set of basis elements generates the subalgebra. Thus you can think of fold() and unfold() as renumbering the basis vectors within a frame.
Eg. {-4,-2,1,3}.fold({-4,-2,1,3}) == {-2,-1,1,2}
{-2,3}.fold({-4,-2,1,3}) == {-1,2}
Operations in matrix_multi<> take place in the folded subalgebra.
The common interface for framed_multi<> and matrix_multi<> is defined via the abstract base template class clifford_algebra<Scalar_T, Index_Set_T, Multivector_T>. This uses what I call pseudo-polymorphism. framed_multi<> and matrix_multi<> do not have a common base class, only a common interface defined via the template. The abstract base template class is used only to define the derived classes and is not used in the library. This keeps overhead to a minimum.
Non-member function code common between the two derived template classes is implemented via the use of template template functions.
For example, the functions sqrt, exp, log, cos, cosh, sin, sinh, tan, tanh, acos, acosh, asin, asinh, atan and atanh use algorithms based on matrix algorithms and Pade approximations, in such a way that the same algorithm is used for both framed_multi<> and matrix_multi<>.
The log and sqrt functions are based on Pade approximation and on matrix algorithms recently published by Cheng, Higham, Kenney and Laub.
The exp function is based on Pade approximation, and the other transcendental functions are based on exp and log. References for these are Abramowitz and Stegun, Gerald and Wheatley, Golub and van Loan. See details below.
GluCat uses the C++ Standard Library and the Boost Library. Make sure you have these installed and working before attempting to use GluCat. Use the Boost download page to obtain the Boost Library.
Make sure that you have installed Boost with the same version of g++ as you will be using for GluCat.
The headers need a compiler which correctly handles template template classes.
GluCat 0.4.0 has not been tested with Boost versions before 1.33.0.
With g++ 4.0.0 and g++ 4.0.1, optimization with -O3 may corrupt iterators. Use -O3 -fno-strict-aliasing instead. See gcc bug 23599.
GluCat 0.4.0 has so far been tested using:
There is plenty of testing, cleanup, packaging and documentation work yet to do.
Transcendental functions still need improvement.
The interface and library could be further improved in a number of different ways.
Paul C. Leopardi
Dr. William McLean, School of Mathematics, UNSW.
Arvind Raja
Arvind Raja's prototype was modelled on Clical, by Pertti Lounesto, R. Mikkola, V. Vierros.
Joerg Walter.
Joerg Walter, John Fletcher, Henk Jansen, Johannes Brunen.
John Fletcher: BoostCliffordDiscussion
Carlos O'Ryan: workaround used in generator.h for: "... only defines a private destructor and has no friends"
| |