In order to implement Lie map methods in MARYLIE, we must efficiently represent and symbolically calculate with truncated power series. Analytically we can represent a monomial in a 2N dimensional phase space as a coefficient a and a vector j of 2N non-negative integer entries giving the exponents of the phase space variables in a fixed order, e.g. for N = 1, we have [3.14,(2,5)] <==> 3.14x2p5. The obvious corresponding representation on the computer is to store coefficients of a polynomial in a table indexed by the entries of the exponent vector of its constituent monomials:
double p[MAXDEG+1][MAXDEG+1];
p[2][5] = 3.14;
While this does work in principle, it quickly becomes clear on inspection that
if we restrict ourselves to polynomials of degree less than some fixed degree,
which we always do in MARYLIE, then many of the allocated locations in the
array p will always be unused. For high accuracy calculations,
this wasted memory becomes prohibitively large.
To circumvent this memory waste problem, we fix the maximum degree d of polynomials in advance, and use a bijective mapping Z+2N <--> Z+ that puts the monomials in a graded lexicographic ordering. This allows us to work with a singly indexed array in which all entries correspond to a monomial m with deg m < d+1.
The Lie methods implemented in MARYLIE require us to perform three basic operations on phase space polynomials: addition, multiplication, and Poisson bracketing. The first, addition, is trivially implemented as vector addition between two polynomial arrays. The remaining two operations are less trivial, and can be quite computationally intensive. Because MARYLIE needs to compute many products and Poisson brackets, it is highly advantageous to optimize these calculations. Our current scheme, in fact, eliminates nearly all runtime calculations associated with performing these operations by storing the results of the required index manipulations for a general multiplication or Poisson bracket in look up tables. At runtime, then, the algorithm makes a call to memory rather than calculating the result, essentially trading storage memory space for computational speed.
The current work is done primarily in FORTRAN77. To verify the correctness of our algorithms, I also code in F90 to utilize the LBNL version of Martin Berz's TPSA package with operator overloading by Etienne Forest and Frank Schmidt, and in C/C++ for other general purpose routines.
For our work, however, we would like the Clebsch-Gordan series for tensor products between representations of the symplectic groups. This would allow us, for example, to project out the Hamiltonian piece of a vector field be Lie algebraic methods. I am currently working on a method for calculating the Clebsch-Gordan series by formulating raising and lowering operators to act on the basis states of general representation labeled in a manner due to Klimyk. Although I am currently restricting my attention to representations of Sp(4), I hope to be able to generalize my results to Sp(2N) by an inductive method.
As a toy problem for the proposed induction, I recently proved that the dimension of a representation specified by a highest weight w = (w1, ... , wN) can be expressed as a polynomial in the dimensions of the Sp(2) representations specified by the one dimensional highest weight vectors (i.e. scalars) wi, i = 1, ..., N. I am currently preparing this result for publication.