| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
By reading so far you should have gotten a fairly good understanding of GiNaC's design patterns. From here on you should start reading the sources. All we can do now is issue some recommendations how to tackle GiNaC's many loose ends in order to fulfill everybody's dreams. If you develop some useful extension please don't hesitate to contact the GiNaC authors--they will happily incorporate them into future versions.
| 6.1 What doesn't belong into GiNaC | What to avoid. | |
| 6.2 Symbolic functions | Implementing symbolic functions. | |
| 6.3 GiNaC's expression output system | Adding new output formats. | |
| 6.4 Structures | Defining new algebraic classes (the easy way). | |
| 6.5 Adding classes | Defining new algebraic classes (the hard way). |
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
First of all, GiNaC's name must be read literally. It is designed to be
a library for use within C++. The tiny ginsh accompanying
GiNaC makes this even more clear: it doesn't even attempt to provide a
language. There are no loops or conditional expressions in
ginsh, it is merely a window into the library for the
programmer to test stuff (or to show off). Still, the design of a
complete CAS with a language of its own, graphical capabilities and all
this on top of GiNaC is possible and is without doubt a nice project for
the future.
There are many built-in functions in GiNaC that do not know how to
evaluate themselves numerically to a precision declared at runtime
(using Digits). Some may be evaluated at certain points, but not
generally. This ought to be fixed. However, doing numerical
computations with GiNaC's quite abstract classes is doomed to be
inefficient. For this purpose, the underlying foundation classes
provided by CLN are much better suited.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
The easiest and most instructive way to start extending GiNaC is probably to create your own symbolic functions. These are implemented with the help of two preprocessor macros:
DECLARE_FUNCTION_<n>P(<name>) REGISTER_FUNCTION(<name>, <options>) |
The DECLARE_FUNCTION macro will usually appear in a header file. It
declares a C++ function with the given `name' that takes exactly `n'
parameters of type ex and returns a newly constructed GiNaC
function object that represents your function.
The REGISTER_FUNCTION macro implements the function. It must be passed
the same `name' as the respective DECLARE_FUNCTION macro, and a
set of options that associate the symbolic function with C++ functions you
provide to implement the various methods such as evaluation, derivative,
series expansion etc. They also describe additional attributes the function
might have, such as symmetry and commutation properties, and a name for
LaTeX output. Multiple options are separated by the member access operator
`.' and can be given in an arbitrary order.
(By the way: in case you are worrying about all the macros above we can assure you that functions are GiNaC's most macro-intense classes. We have done our best to avoid macros where we can.)
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
Here is an example for the implementation of a function with two arguments that is not further evaluated:
DECLARE_FUNCTION_2P(myfcn) REGISTER_FUNCTION(myfcn, dummy()) |
Any code that has seen the DECLARE_FUNCTION line can use myfcn()
in algebraic expressions:
{
...
symbol x("x");
ex e = 2*myfcn(42, 1+3*x) - x;
cout << e << endl;
// prints '2*myfcn(42,1+3*x)-x'
...
}
|
The dummy() option in the REGISTER_FUNCTION line signifies
"no options". A function with no options specified merely acts as a kind of
container for its arguments. It is a pure "dummy" function with no associated
logic (which is, however, sometimes perfectly sufficient).
Let's now have a look at the implementation of GiNaC's cosine function for an example of how to make an "intelligent" function.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
The GiNaC header file `inifcns.h' contains the line
DECLARE_FUNCTION_1P(cos) |
which declares to all programs using GiNaC that there is a function `cos'
that takes one ex as an argument. This is all they need to know to use
this function in expressions.
The implementation of the cosine function is in `inifcns_trans.cpp'. Here
is its REGISTER_FUNCTION line:
REGISTER_FUNCTION(cos, eval_func(cos_eval).
evalf_func(cos_evalf).
derivative_func(cos_deriv).
latex_name("\\cos"));
|
There are four options defined for the cosine function. One of them
(latex_name) gives the function a proper name for LaTeX output; the
other three indicate the C++ functions in which the "brains" of the cosine
function are defined.
The eval_func() option specifies the C++ function that implements
the eval() method, GiNaC's anonymous evaluator. This function takes
the same number of arguments as the associated symbolic function (one in this
case) and returns the (possibly transformed or in some way simplified)
symbolically evaluated function (See section 4.2 Automatic evaluation and canonicalization of expressions, for a description
of the automatic evaluation process). If no (further) evaluation is to take
place, the eval_func() function must return the original function
with .hold(), to avoid a potential infinite recursion. If your
symbolic functions produce a segmentation fault or stack overflow when
using them in expressions, you are probably missing a .hold()
somewhere.
The eval_func() function for the cosine looks something like this
(actually, it doesn't look like this at all, but it should give you an idea
what is going on):
static ex cos_eval(const ex & x)
{
if ("x is a multiple of 2*Pi")
return 1;
else if ("x is a multiple of Pi")
return -1;
else if ("x is a multiple of Pi/2")
return 0;
// more rules...
else if ("x has the form 'acos(y)'")
return y;
else if ("x has the form 'asin(y)'")
return sqrt(1-y^2);
// more rules...
else
return cos(x).hold();
}
|
This function is called every time the cosine is used in a symbolic expression:
{
...
e = cos(Pi);
// this calls cos_eval(Pi), and inserts its return value into
// the actual expression
cout << e << endl;
// prints '-1'
...
}
|
In this way, cos(4*Pi) automatically becomes ,
cos(asin(a+b)) becomes sqrt(1-(a+b)^2), etc. If no reasonable
symbolic transformation can be done, the unmodified function is returned
with .hold().
GiNaC doesn't automatically transform cos(2) to `-0.416146...'.
The user has to call evalf() for that. This is implemented in a
different function:
static ex cos_evalf(const ex & x)
{
if (is_a<numeric>(x))
return cos(ex_to<numeric>(x));
else
return cos(x).hold();
}
|
Since we are lazy we defer the problem of numeric evaluation to somebody else,
in this case the cos() function for numeric objects, which in
turn hands it over to the cos() function in CLN. The .hold()
isn't really needed here, but reminds us that the corresponding eval()
function would require it in this place.
Differentiation will surely turn up and so we need to tell cos
what its first derivative is (higher derivatives, .diff(x,3) for
instance, are then handled automatically by basic::diff and
ex::diff):
static ex cos_deriv(const ex & x, unsigned diff_param)
{
return -sin(x);
}
|
The second parameter is obligatory but uninteresting at this point. It specifies which parameter to differentiate in a partial derivative in case the function has more than one parameter, and its main application is for correct handling of the chain rule.
An implementation of the series expansion is not needed for cos() as
it doesn't have any poles and GiNaC can do Taylor expansion by itself (as
long as it knows what the derivative of cos() is). tan(), on
the other hand, does have poles and may need to do Laurent expansion:
static ex tan_series(const ex & x, const relational & rel,
int order, unsigned options)
{
// Find the actual expansion point
const ex x_pt = x.subs(rel);
if ("x_pt is not an odd multiple of Pi/2")
throw do_taylor(); // tell function::series() to do Taylor expansion
// On a pole, expand sin()/cos()
return (sin(x)/cos(x)).series(rel, order+2, options);
}
|
The series() implementation of a function must return a
pseries object, otherwise your code will crash.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
GiNaC functions understand several more options which are always
specified as .option(params). None of them are required, but you
need to specify at least one option to REGISTER_FUNCTION(). There
is a do-nothing option called dummy() which you can use to define
functions without any special options.
eval_func(<C++ function>) evalf_func(<C++ function>) derivative_func(<C++ function>) series_func(<C++ function>) conjugate_func(<C++ function>) |
These specify the C++ functions that implement symbolic evaluation,
numeric evaluation, partial derivatives, and series expansion, respectively.
They correspond to the GiNaC methods eval(), evalf(),
diff() and series().
The eval_func() function needs to use .hold() if no further
automatic evaluation is desired or possible.
If no series_func() is given, GiNaC defaults to simple Taylor
expansion, which is correct if there are no poles involved. If the function
has poles in the complex plane, the series_func() needs to check
whether the expansion point is on a pole and fall back to Taylor expansion
if it isn't. Otherwise, the pole usually needs to be regularized by some
suitable transformation.
latex_name(const string & n) |
specifies the LaTeX code that represents the name of the function in LaTeX
output. The default is to put the function name in an \mbox{}.
do_not_evalf_params() |
This tells evalf() to not recursively evaluate the parameters of the
function before calling the evalf_func().
set_return_type(unsigned return_type, unsigned return_type_tinfo) |
This allows you to explicitly specify the commutation properties of the
function (See section 4.15 Non-commutative objects, for an explanation of
(non)commutativity in GiNaC). For example, you can use
set_return_type(return_types::noncommutative, TINFO_matrix) to make
GiNaC treat your function like a matrix. By default, functions inherit the
commutation properties of their first argument.
set_symmetry(const symmetry & s) |
specifies the symmetry properties of the function with respect to its arguments. See section 4.14 Indexed objects, for an explanation of symmetry specifications. GiNaC will automatically rearrange the arguments of symmetric functions into a canonical order.
Sometimes you may want to have finer control over how functions are
displayed in the output. For example, the abs() function prints
itself as `abs(x)' in the default output format, but as `|x|'
in LaTeX mode, and fabs(x) in C source output. This is achieved
with the
print_func<C>(<C++ function>) |
option which is explained in the next section.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
The DECLARE_FUNCTION and REGISTER_FUNCTION macros define
functions with a fixed number of arguments. Sometimes, though, you may need
to have a function that accepts a variable number of expressions. One way to
accomplish this is to pass variable-length lists as arguments. The
Li() function uses this method for multiple polylogarithms.
It is also possible to define functions that accept a different number of
parameters under the same function name, such as the psi() function
which can be called either as psi(z) (the digamma function) or as
psi(n, z) (polygamma functions). These are actually two different
functions in GiNaC that, however, have the same name. Defining such
functions is not possible with the macros but requires manually fiddling
with GiNaC internals. If you are interested, please consult the GiNaC source
code for the psi() function (`inifcns.h' and
`inifcns_gamma.cpp').
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
GiNaC allows the output of expressions in a variety of different formats (see section 5.15 Input and output of expressions). This section will explain how expression output is implemented internally, and how to define your own output formats or change the output format of built-in algebraic objects. You will also want to read this section if you plan to write your own algebraic classes or functions.
All the different output formats are represented by a hierarchy of classes
rooted in the print_context class, defined in the `print.h'
header file:
print_dflt
print_latex
print_tree
print_csrc
print_csrc_float
float type
print_csrc_double
double type
print_csrc_cl_N
The print_context base class provides two public data members:
class print_context
{
...
public:
std::ostream & s;
unsigned options;
};
|
s is a reference to the stream to output to, while options
holds flags and modifiers. Currently, there is only one flag defined:
print_options::print_index_dimensions instructs the idx class
to print the index dimension which is normally hidden.
When you write something like std::cout << e, where e is
an object of class ex, GiNaC will construct an appropriate
print_context object (of a class depending on the selected output
format), fill in the s and options members, and call
void ex::print(const print_context & c, unsigned level = 0) const; |
which in turn forwards the call to the print() method of the
top-level algebraic object contained in the expression.
Unlike other methods, GiNaC classes don't usually override their
print() method to implement expression output. Instead, the default
implementation basic::print(c, level) performs a run-time double
dispatch to a function selected by the dynamic type of the object and the
passed print_context. To this end, GiNaC maintains a separate method
table for each class, similar to the virtual function table used for ordinary
(single) virtual function dispatch.
The method table contains one slot for each possible print_context
type, indexed by the (internally assigned) serial number of the type. Slots
may be empty, in which case GiNaC will retry the method lookup with the
print_context object's parent class, possibly repeating the process
until it reaches the print_context base class. If there's still no
method defined, the method table of the algebraic object's parent class
is consulted, and so on, until a matching method is found (eventually it
will reach the combination basic/print_context, which prints the
object's class name enclosed in square brackets).
You can think of the print methods of all the different classes and output
formats as being arranged in a two-dimensional matrix with one axis listing
the algebraic classes and the other axis listing the print_context
classes.
Subclasses of basic can, of course, also overload basic::print()
to implement printing, but then they won't get any of the benefits of the
double dispatch mechanism (such as the ability for derived classes to
inherit only certain print methods from its parent, or the replacement of
methods at run-time).
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
The method table for a class is set up either in the definition of the class,
by passing the appropriate print_func<C>() option to
GINAC_IMPLEMENT_REGISTERED_CLASS_OPT() (See section 6.5 Adding classes, for
an example), or at run-time using set_print_func<T, C>(). The latter
can also be used to override existing methods dynamically.
The argument to print_func<C>() and set_print_func<T, C>() can
be a member function of the class (or one of its parent classes), a static
member function, or an ordinary (global) C++ function. The C template
parameter specifies the appropriate print_context type for which the
method should be invoked, while, in the case of set_print_func<>(), the
T parameter specifies the algebraic class (for print_func<>(),
the class is the one being implemented by
GINAC_IMPLEMENT_REGISTERED_CLASS_OPT).
For print methods that are member functions, their first argument must be of
a type convertible to a const C &, and the second argument must be an
unsigned.
For static members and global functions, the first argument must be of a type
convertible to a const T &, the second argument must be of a type
convertible to a const C &, and the third argument must be an
unsigned. A global function will, of course, not have access to
private and protected members of T.
The unsigned argument of the print methods (and of ex::print()
and basic::print()) is used for proper parenthesizing of the output
(and by print_tree for proper indentation). It can be used for similar
purposes if you write your own output formats.
The explanations given above may seem complicated, but in practice it's really simple, as shown in the following example. Suppose that we want to display exponents in LaTeX output not as superscripts but with little upwards-pointing arrows. This can be achieved in the following way:
void my_print_power_as_latex(const power & p,
const print_latex & c,
unsigned level)
{
// get the precedence of the 'power' class
unsigned power_prec = p.precedence();
// if the parent operator has the same or a higher precedence
// we need parentheses around the power
if (level >= power_prec)
c.s << '(';
// print the basis and exponent, each enclosed in braces, and
// separated by an uparrow
c.s << '{';
p.op(0).print(c, power_prec);
c.s << "}\\uparrow{";
p.op(1).print(c, power_prec);
c.s << '}';
// don't forget the closing parenthesis
if (level >= power_prec)
c.s << ')';
}
int main()
{
// a sample expression
symbol x("x"), y("y");
ex e = -3*pow(x, 3)*pow(y, -2) + pow(x+y, 2) - 1;
// switch to LaTeX mode
cout << latex;
// this prints "-1+{(y+x)}^{2}-3 \frac{x^{3}}{y^{2}}"
cout << e << endl;
// now we replace the method for the LaTeX output of powers with
// our own one
set_print_func<power, print_latex>(my_print_power_as_latex);
// this prints "-1+{{(y+x)}}\uparrow{2}-3 \frac{{x}\uparrow{3}}{{y}
// \uparrow{2}}"
cout << e << endl;
}
|
Some notes:
my_print_power_as_latex could also have been
a const basic &, the second one a const print_context &.
mul objects converting their operands to
power objects for the purpose of printing.
mul class.
power/print_latex method provided by GiNaC prints square roots
using \sqrt, but the above code doesn't.
It's not possible to restore a method table entry to its previous or default
value. Once you have called set_print_func(), you can only override
it with another call to set_print_func(), but you can't easily go back
to the default behavior again (you can, of course, dig around in the GiNaC
sources, find the method that is installed at startup
(power::do_print_latex in this case), and set_print_func that
one; that is, after you circumvent the C++ member access control...).
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
Symbolic functions employ a print method dispatch mechanism similar to the
one used for classes. The methods are specified with print_func<C>()
function options. If you don't specify any special print methods, the function
will be printed with its name (or LaTeX name, if supplied), followed by a
comma-separated list of arguments enclosed in parentheses.
For example, this is what GiNaC's `abs()' function is defined like:
static ex abs_eval(const ex & arg) { ... }
static ex abs_evalf(const ex & arg) { ... }
static void abs_print_latex(const ex & arg, const print_context & c)
{
c.s << "{|"; arg.print(c); c.s << "|}";
}
static void abs_print_csrc_float(const ex & arg, const print_context & c)
{
c.s << "fabs("; arg.print(c); c.s << ")";
}
REGISTER_FUNCTION(abs, eval_func(abs_eval).
evalf_func(abs_evalf).
print_func<print_latex>(abs_print_latex).
print_func<print_csrc_float>(abs_print_csrc_float).
print_func<print_csrc_double>(abs_print_csrc_float));
|
This will display `abs(x)' as `|x|' in LaTeX mode and fabs(x)
in non-CLN C source output, but as abs(x) in all other formats.
There is currently no equivalent of set_print_func() for functions.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
Creating a new output format involves subclassing print_context,
which is somewhat similar to adding a new algebraic class
(see section 6.5 Adding classes). There is a macro GINAC_DECLARE_PRINT_CONTEXT
that needs to go into the class definition, and a corresponding macro
GINAC_IMPLEMENT_PRINT_CONTEXT that has to appear at global scope.
Every print_context class needs to provide a default constructor
and a constructor from an std::ostream and an unsigned
options value.
Here is an example for a user-defined print_context class:
class print_myformat : public print_dflt
{
GINAC_DECLARE_PRINT_CONTEXT(print_myformat, print_dflt)
public:
print_myformat(std::ostream & os, unsigned opt = 0)
: print_dflt(os, opt) {}
};
print_myformat::print_myformat() : print_dflt(std::cout) {}
GINAC_IMPLEMENT_PRINT_CONTEXT(print_myformat, print_dflt)
|
That's all there is to it. None of the actual expression output logic is implemented in this class. It merely serves as a selector for choosing a particular format. The algorithms for printing expressions in the new format are implemented as print methods, as described above.
print_myformat is a subclass of print_dflt, so it behaves
exactly like GiNaC's default output format:
{
symbol x("x");
ex e = pow(x, 2) + 1;
// this prints "1+x^2"
cout << e << endl;
// this also prints "1+x^2"
e.print(print_myformat()); cout << endl;
...
}
|
To fill print_myformat with life, we need to supply appropriate
print methods with set_print_func(), like this:
// This prints powers with '**' instead of '^'. See the LaTeX output
// example above for explanations.
void print_power_as_myformat(const power & p,
const print_myformat & c,
unsigned level)
{
unsigned power_prec = p.precedence();
if (level >= power_prec)
c.s << '(';
p.op(0).print(c, power_prec);
c.s << "**";
p.op(1).print(c, power_prec);
if (level >= power_prec)
c.s << ')';
}
{
...
// install a new print method for power objects
set_print_func<power, print_myformat>(print_power_as_myformat);
// now this prints "1+x**2"
e.print(print_myformat()); cout << endl;
// but the default format is still "1+x^2"
cout << e << endl;
}
|
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
If you are doing some very specialized things with GiNaC, or if you just
need some more organized way to store data in your expressions instead of
anonymous lists, you may want to implement your own algebraic classes.
('algebraic class' means any class directly or indirectly derived from
basic that can be used in GiNaC expressions).
GiNaC offers two ways of accomplishing this: either by using the
structure<T> template class, or by rolling your own class from
scratch. This section will discuss the structure<T> template which
is easier to use but more limited, while the implementation of custom
GiNaC classes is the topic of the next section. However, you may want to
read both sections because many common concepts and member functions are
shared by both concepts, and it will also allow you to decide which approach
is most suited to your needs.
The structure<T> template, defined in the GiNaC header file
`structure.h', wraps a type that you supply (usually a C++ struct
or class) into a GiNaC object that can be used in expressions.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
Let's suppose that we need a way to handle some kind of abstract scalar
product of the form `<x|y>' in expressions. Objects of the scalar
product class have to store their left and right operands, which can in turn
be arbitrary expressions. Here is a possible way to represent such a
product in a C++ struct:
#include <iostream>
using namespace std;
#include <ginac/ginac.h>
using namespace GiNaC;
struct sprod_s {
ex left, right;
sprod_s() {}
sprod_s(ex l, ex r) : left(l), right(r) {}
};
|
The default constructor is required. Now, to make a GiNaC class out of this data structure, we need only one line:
typedef structure<sprod_s> sprod; |
That's it. This line constructs an algebraic class sprod which
contains objects of type sprod_s. We can now use sprod in
expressions like any other GiNaC class:
...
symbol a("a"), b("b");
ex e = sprod(sprod_s(a, b));
...
|
Note the difference between sprod which is the algebraic class, and
sprod_s which is the unadorned C++ structure containing the left
and right data members. As shown above, an sprod can be
constructed from an sprod_s object.
If you find the nested sprod(sprod_s()) constructor too unwieldy,
you could define a little wrapper function like this:
inline ex make_sprod(ex left, ex right)
{
return sprod(sprod_s(left, right));
}
|
The sprod_s object contained in sprod can be accessed with
the GiNaC ex_to<>() function followed by the -> operator or
get_struct():
...
cout << ex_to<sprod>(e)->left << endl;
// -> a
cout << ex_to<sprod>(e).get_struct().right << endl;
// -> b
...
|
You only have read access to the members of sprod_s.
The type definition of sprod is enough to write your own algorithms
that deal with scalar products, for example:
ex swap_sprod(ex p)
{
if (is_a<sprod>(p)) {
const sprod_s & sp = ex_to<sprod>(p).get_struct();
return make_sprod(sp.right, sp.left);
} else
return p;
}
...
f = swap_sprod(e);
// f is now <b|a>
...
|
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
While the sprod type is useable it still leaves something to be
desired, most notably proper output:
...
cout << e << endl;
// -> [structure object]
...
|
By default, any structure types you define will be printed as
`[structure object]'. To override this you can either specialize the
template's print() member function, or specify print methods with
set_print_func<>(), as described in 6.3 GiNaC's expression output system. Unfortunately,
it's not possible to supply class options like print_func<>() to
structures, so for a self-contained structure type you need to resort to
overriding the print() function, which is also what we will do here.
The member functions of GiNaC classes are described in more detail in the next section, but it shouldn't be hard to figure out what's going on here:
void sprod::print(const print_context & c, unsigned level) const
{
// tree debug output handled by superclass
if (is_a<print_tree>(c))
inherited::print(c, level);
// get the contained sprod_s object
const sprod_s & sp = get_struct();
// print_context::s is a reference to an ostream
c.s << "<" << sp.left << "|" << sp.right << ">";
}
|
Now we can print expressions containing scalar products:
...
cout << e << endl;
// -> <a|b>
cout << swap_sprod(e) << endl;
// -> <b|a>
...
|
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
The sprod class defined so far still has one important drawback: all
scalar products are treated as being equal because GiNaC doesn't know how to
compare objects of type sprod_s. This can lead to some confusing
and undesired behavior:
...
cout << make_sprod(a, b) - make_sprod(a*a, b*b) << endl;
// -> 0
cout << make_sprod(a, b) + make_sprod(a*a, b*b) << endl;
// -> 2*<a|b> or 2*<a^2|b^2> (which one is undefined)
...
|
To remedy this, we first need to define the operators == and <
for objects of type sprod_s:
inline bool operator==(const sprod_s & lhs, const sprod_s & rhs)
{
return lhs.left.is_equal(rhs.left) && lhs.right.is_equal(rhs.right);
}
inline bool operator<(const sprod_s & lhs, const sprod_s & rhs)
{
return lhs.left.compare(rhs.left) < 0
? true : lhs.right.compare(rhs.right) < 0;
}
|
The ordering established by the < operator doesn't have to make any
algebraic sense, but it needs to be well defined. Note that we can't use
expressions like lhs.left == rhs.left or lhs.left < rhs.left
in the implementation of these operators because they would construct
GiNaC relational objects which in the case of < do not
establish a well defined ordering (for arbitrary expressions, GiNaC can't
decide which one is algebraically 'less').
Next, we need to change our definition of the sprod type to let
GiNaC know that an ordering relation exists for the embedded objects:
typedef structure<sprod_s, compare_std_less> sprod; |
sprod objects then behave as expected:
...
cout << make_sprod(a, b) - make_sprod(a*a, b*b) << endl;
// -> <a|b>-<a^2|b^2>
cout << make_sprod(a, b) + make_sprod(a*a, b*b) << endl;
// -> <a|b>+<a^2|b^2>
cout << make_sprod(a, b) - make_sprod(a, b) << endl;
// -> 0
cout << make_sprod(a, b) + make_sprod(a, b) << endl;
// -> 2*<a|b>
...
|
The compare_std_less policy parameter tells GiNaC to use the
std::less and std::equal_to functors to compare objects of
type sprod_s. By default, these functors forward their work to the
standard < and == operators, which we have overloaded.
Alternatively, we could have specialized std::less and
std::equal_to for class sprod_s.
GiNaC provides two other comparison policies for structure<T>
objects: the default compare_all_equal, and compare_bitwise
which does a bit-wise comparison of the contained T objects.
This should be used with extreme care because it only works reliably with
built-in integral types, and it also compares any padding (filler bytes of
undefined value) that the T class might have.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
Our scalar product class has two subexpressions: the left and right
operands. It might be a good idea to make them accessible via the standard
nops() and op() methods:
size_t sprod::nops() const
{
return 2;
}
ex sprod::op(size_t i) const
{
switch (i) {
case 0:
return get_struct().left;
case 1:
return get_struct().right;
default:
throw std::range_error("sprod::op(): no such operand");
}
}
|
Implementing nops() and op() for container types such as
sprod has two other nice side effects:
has() works as expected
calchash() takes subexpressions into account)
There is a non-const variant of op() called let_op() that
allows replacing subexpressions:
ex & sprod::let_op(size_t i)
{
// every non-const member function must call this
ensure_if_modifiable();
switch (i) {
case 0:
return get_struct().left;
case 1:
return get_struct().right;
default:
throw std::range_error("sprod::let_op(): no such operand");
}
}
|
Once we have provided let_op() we also get subs() and
map() for free. In fact, every container class that returns a non-null
nops() value must either implement let_op() or provide custom
implementations of subs() and map().
In turn, the availability of map() enables the recursive behavior of a
couple of other default method implementations, in particular evalf(),
evalm(), normal(), diff() and expand(). Although
we probably want to provide our own version of expand() for scalar
products that turns expressions like `<a+b|c>' into `<a|c>+<b|c>'.
This is left as an exercise for the reader.
The structure<T> template defines many more member functions that
you can override by specialization to customize the behavior of your
structures. You are referred to the next section for a description of
some of these (especially eval()). There is, however, one topic
that shall be addressed here, as it demonstrates one peculiarity of the
structure<T> template: archiving.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
If you don't know how the archiving of GiNaC objects is implemented, you should first read the next section and then come back here. You're back? Good.
To implement archiving for structures it is not enough to provide
specializations for the archive() member function and the
unarchiving constructor (the unarchive() function has a default
implementation). You also need to provide a unique name (as a string literal)
for each structure type you define. This is because in GiNaC archives,
the class of an object is stored as a string, the class name.
By default, this class name (as returned by the class_name() member
function) is `structure' for all structure classes. This works as long
as you have only defined one structure type, but if you use two or more you
need to provide a different name for each by specializing the
get_class_name() member function. Here is a sample implementation
for enabling archiving of the scalar product type defined above:
const char *sprod::get_class_name() { return "sprod"; }
void sprod::archive(archive_node & n) const
{
inherited::archive(n);
n.add_ex("left", get_struct().left);
n.add_ex("right", get_struct().right);
}
sprod::structure(const archive_node & n, lst & sym_lst) : inherited(n, sym_lst)
{
n.find_ex("left", get_struct().left, sym_lst);
n.find_ex("right", get_struct().right, sym_lst);
}
|
Note that the unarchiving constructor is sprod::structure and not
sprod::sprod, and that we don't need to supply an
sprod::unarchive() function.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
The structure<T> template provides an way to extend GiNaC with custom
algebraic classes that is easy to use but has its limitations, the most
severe of which being that you can't add any new member functions to
structures. To be able to do this, you need to write a new class definition
from scratch.
This section will explain how to implement new algebraic classes in GiNaC by giving the example of a simple 'string' class. After reading this section you will know how to properly declare a GiNaC class and what the minimum required member functions are that you have to implement. We only cover the implementation of a 'leaf' class here (i.e. one that doesn't contain subexpressions). Creating a container class like, for example, a class representing tensor products is more involved but this section should give you enough information so you can consult the source to GiNaC's predefined classes if you want to implement something more complicated.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
All algebraic classes (that is, all classes that can appear in expressions)
in GiNaC are direct or indirect subclasses of the class basic. So a
basic * (which is essentially what an ex is) represents a
generic pointer to an algebraic class. Occasionally it is necessary to find
out what the class of an object pointed to by a basic * really is.
Also, for the unarchiving of expressions it must be possible to find the
unarchive() function of a class given the class name (as a string). A
system that provides this kind of information is called a run-time type
information (RTTI) system. The C++ language provides such a thing (see the
standard header file `<typeinfo>') but for efficiency reasons GiNaC
implements its own, simpler RTTI.
The RTTI in GiNaC is based on two mechanisms:
basic class declares a member variable tinfo_key which
holds an unsigned integer that identifies the object's class. These numbers
are defined in the `tinfos.h' header file for the built-in GiNaC
classes. They all start with TINFO_.
basic. The information
available includes the class names, the tinfo_keys, and pointers
to the unarchiving functions. This class registry is defined in the
`registrar.h' header file.
The disadvantage of this proprietary RTTI implementation is that there's a little more to do when implementing new classes (C++'s RTTI works more or less automatically) but don't worry, most of the work is simplified by macros.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
Now we will start implementing a new class mystring that allows
placing character strings in algebraic expressions (this is not very useful,
but it's just an example). This class will be a direct subclass of
basic. You can use this sample implementation as a starting point
for your own classes.
The code snippets given here assume that you have included some header files as follows:
#include <iostream> #include <string> #include <stdexcept> using namespace std; #include <ginac/ginac.h> using namespace GiNaC; |
The first thing we have to do is to define a tinfo_key for our new
class. This can be any arbitrary unsigned number that is not already taken
by one of the existing classes but it's better to come up with something
that is unlikely to clash with keys that might be added in the future. The
numbers in `tinfos.h' are modeled somewhat after the class hierarchy
which is not a requirement but we are going to stick with this scheme:
const unsigned TINFO_mystring = 0x42420001U; |
Now we can write down the class declaration. The class stores a C++
string and the user shall be able to construct a mystring
object from a C or C++ string:
class mystring : public basic
{
GINAC_DECLARE_REGISTERED_CLASS(mystring, basic)
public:
mystring(const string &s);
mystring(const char *s);
private:
string str;
};
GINAC_IMPLEMENT_REGISTERED_CLASS(mystring, basic)
|
The GINAC_DECLARE_REGISTERED_CLASS and GINAC_IMPLEMENT_REGISTERED_CLASS
macros are defined in `registrar.h'. They take the name of the class
and its direct superclass as arguments and insert all required declarations
for the RTTI system. The GINAC_DECLARE_REGISTERED_CLASS should be
the first line after the opening brace of the class definition. The
GINAC_IMPLEMENT_REGISTERED_CLASS may appear anywhere else in the
source (at global scope, of course, not inside a function).
GINAC_DECLARE_REGISTERED_CLASS contains, among other things the
declarations of the default constructor and a couple of other functions that
are required. It also defines a type inherited which refers to the
superclass so you don't have to modify your code every time you shuffle around
the class hierarchy. GINAC_IMPLEMENT_REGISTERED_CLASS registers the
class with the GiNaC RTTI (there is also a
GINAC_IMPLEMENT_REGISTERED_CLASS_OPT which allows specifying additional
options for the class, and which we will be using instead in a few minutes).
Now there are seven member functions we have to implement to get a working class:
mystring(), the default constructor.
void archive(archive_node &n), the archiving function. This stores all
information needed to reconstruct an object of this class inside an
archive_node.
mystring(const archive_node &n, lst &sym_lst), the unarchiving
constructor. This constructs an instance of the class from the information
found in an archive_node.
ex unarchive(const archive_node &n, lst &sym_lst), the static
unarchiving function. It constructs a new instance by calling the unarchiving
constructor.
int compare_same_type(const basic &other), which is used internally
by GiNaC to establish a canonical sort order for terms. It returns 0, +1 or
-1, depending on the relative order of this object and the other
object. If it returns 0, the objects are considered equal.
Please notice: This has nothing to do with the (numeric) ordering
relationship expressed by <, >= etc (which cannot be defined
for non-numeric classes). For example, numeric(1).compare_same_type(numeric(2))
may return +1 even though 1 is clearly smaller than 2. Every GiNaC class
must provide a compare_same_type() function, even those representing
objects for which no reasonable algebraic ordering relationship can be
defined.
mystring(const string &s) and mystring(const char *s)
which are the two constructors we declared.
Let's proceed step-by-step. The default constructor looks like this:
mystring::mystring() : inherited(TINFO_mystring) {}
|
The golden rule is that in all constructors you have to set the
tinfo_key member to the TINFO_* value of your class. Otherwise
it will be set by the constructor of the superclass and all hell will break
loose in the RTTI. For your convenience, the basic class provides
a constructor that takes a tinfo_key value, which we are using here
(remember that in our case inherited == basic). If the superclass
didn't have such a constructor, we would have to set the tinfo_key
to the right value manually.
In the default constructor you should set all other member variables to
reasonable default values (we don't need that here since our str
member gets set to an empty string automatically).
Next are the three functions for archiving. You have to implement them even if you don't plan to use archives, but the minimum required implementation is really simple. First, the archiving function:
void mystring::archive(archive_node &n) const
{
inherited::archive(n);
n.add_string("string", str);
}
|
The only thing that is really required is calling the archive()
function of the superclass. Optionally, you can store all information you
deem necessary for representing the object into the passed
archive_node. We are just storing our string here. For more
information on how the archiving works, consult the `archive.h' header
file.
The unarchiving constructor is basically the inverse of the archiving function:
mystring::mystring(const archive_node &n, lst &sym_lst) : inherited(n, sym_lst)
{
n.find_string("string", str);
}
|
If you don't need archiving, just leave this function empty (but you must
invoke the unarchiving constructor of the superclass). Note that we don't
have to set the tinfo_key here because it is done automatically
by the unarchiving constructor of the basic class.
Finally, the unarchiving function:
ex mystring::unarchive(const archive_node &n, lst &sym_lst)
{
return (new mystring(n, sym_lst))->setflag(status_flags::dynallocated);
}
|
You don't have to understand how exactly this works. Just copy these
four lines into your code literally (replacing the class name, of
course). It calls the unarchiving constructor of the class and unless
you are doing something very special (like matching archive_nodes
to global objects) you don't need a different implementation. For those
who are interested: setting the dynallocated flag puts the object
under the control of GiNaC's garbage collection. It will get deleted
automatically once it is no longer referenced.
Our compare_same_type() function uses a provided function to compare
the string members:
int mystring::compare_same_type(const basic &other) const
{
const mystring &o = static_cast<const mystring &>(other);
int cmpval = str.compare(o.str);
if (cmpval == 0)
return 0;
else if (cmpval < 0)
return -1;
else
return 1;
}
|
Although this function takes a basic &, it will always be a reference
to an object of exactly the same class (objects of different classes are not
comparable), so the cast is safe. If this function returns 0, the two objects
are considered equal (in the sense that ), so you should compare
all relevant member variables.
Now the only thing missing is our two new constructors:
mystring::mystring(const string &s) : inherited(TINFO_mystring), str(s) {}
mystring::mystring(const char *s) : inherited(TINFO_mystring), str(s) {}
|
No surprises here. We set the str member from the argument and
remember to pass the right tinfo_key to the basic constructor.
That's it! We now have a minimal working GiNaC class that can store strings in algebraic expressions. Let's confirm that the RTTI works:
ex e = mystring("Hello, world!");
cout << is_a<mystring>(e) << endl;
// -> 1 (true)
cout << e.bp->class_name() << endl;
// -> mystring
|
Obviously it does. Let's see what the expression e looks like:
cout << e << endl; // -> [mystring object] |
Hm, not exactly what we expect, but of course the mystring class
doesn't yet know how to print itself. This can be done either by implementing
the print() member function, or, preferably, by specifying a
print_func<>() class option. Let's say that we want to print the string
surrounded by double quotes:
class mystring : public basic
{
...
protected:
void do_print(const print_context &c, unsigned level = 0) const;
...
};
void mystring::do_print(const print_context &c, unsigned level) const
{
// print_context::s is a reference to an ostream
c.s << '\"' << str << '\"';
}
|
The level argument is only required for container classes to
correctly parenthesize the output.
Now we need to tell GiNaC that mystring objects should use the
do_print() member function for printing themselves. For this, we
replace the line
GINAC_IMPLEMENT_REGISTERED_CLASS(mystring, basic) |
with
GINAC_IMPLEMENT_REGISTERED_CLASS_OPT(mystring, basic, print_func<print_context>(&mystring::do_print)) |
Let's try again to print the expression:
cout << e << endl; // -> "Hello, world!" |
Much better. If we wanted to have mystring objects displayed in a
different way depending on the output format (default, LaTeX, etc.), we
would have supplied multiple print_func<>() options with different
template parameters (print_dflt, print_latex, etc.),
separated by dots. This is similar to the way options are specified for
symbolic functions. See section 6.3 GiNaC's expression output system, for a more in-depth description of the
way expression output is implemented in GiNaC.
The mystring class can be used in arbitrary expressions:
e += mystring("GiNaC rulez");
cout << e << endl;
// -> "GiNaC rulez"+"Hello, world!"
|
(GiNaC's automatic term reordering is in effect here), or even
e = pow(mystring("One string"), 2*sin(Pi-mystring("Another string")));
cout << e << endl;
// -> "One string"^(2*sin(-"Another string"+Pi))
|
Whether this makes sense is debatable but remember that this is only an example. At least it allows you to implement your own symbolic algorithms for your objects.
Note that GiNaC's algebraic rules remain unchanged:
e = mystring("Wow") * mystring("Wow");
cout << e << endl;
// -> "Wow"^2
e = pow(mystring("First")-mystring("Second"), 2);
cout << e.expand() << endl;
// -> -2*"First"*"Second"+"First"^2+"Second"^2
|
There's no way to, for example, make GiNaC's add class perform string
concatenation. You would have to implement this yourself.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
When dealing with objects that are just a little more complicated than the
simple string objects we have implemented, chances are that you will want to
have some automatic simplifications or canonicalizations performed on them.
This is done in the evaluation member function eval(). Let's say that
we wanted all strings automatically converted to lowercase with
non-alphabetic characters stripped, and empty strings removed:
class mystring : public basic
{
...
public:
ex eval(int level = 0) const;
...
};
ex mystring::eval(int level) const
{
string new_str;
for (int i=0; i<str.length(); i++) {
char c = str[i];
if (c >= 'A' && c <= 'Z')
new_str += tolower(c);
else if (c >= 'a' && c <= 'z')
new_str += c;
}
if (new_str.length() == 0)
return 0;
else
return mystring(new_str).hold();
}
|
The level argument is used to limit the recursion depth of the
evaluation. We don't have any subexpressions in the mystring
class so we are not concerned with this. If we had, we would call the
eval() functions of the subexpressions with level - 1 as
the argument if level != 1. The hold() member function
sets a flag in the object that prevents further evaluation. Otherwise
we might end up in an endless loop. When you want to return the object
unmodified, use return this->hold();.
Let's confirm that it works:
ex e = mystring("Hello, world!") + mystring("!?#");
cout << e << endl;
// -> "helloworld"
e = mystring("Wow!") + mystring("WOW") + mystring(" W ** o ** W");
cout << e << endl;
// -> 3*"wow"
|
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
We have implemented only a small set of member functions to make the class work in the GiNaC framework. There are two functions that are not strictly required but will make operations with objects of the class more efficient:
unsigned calchash() const; bool is_equal_same_type(const basic &other) const; |
The calchash() method returns an unsigned hash value for the
object which will allow GiNaC to compare and canonicalize expressions much
more efficiently. You should consult the implementation of some of the built-in
GiNaC classes for examples of hash functions. The default implementation of
calchash() calculates a hash value out of the tinfo_key of the
class and all subexpressions that are accessible via op().
is_equal_same_type() works like compare_same_type() but only
tests for equality without establishing an ordering relation, which is often
faster. The default implementation of is_equal_same_type() just calls
compare_same_type() and tests its result for zero.
| [ < ] | [ > ] | [ << ] | [ Up ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
For a real algebraic class, there are probably some more functions that you might want to provide:
bool info(unsigned inf) const; ex evalf(int level = 0) const; ex series(const relational & r, int order, unsigned options = 0) const; ex derivative(const symbol & s) const; |
If your class stores sub-expressions (see the scalar product example in the previous section) you will probably want to override
size_t nops() cont; ex op(size_t i) const; ex & let_op(size_t i); ex subs(const lst & ls, const lst & lr, unsigned options = 0) const; ex map(map_function & f) const; |
let_op() is a variant of op() that allows write access. The
default implementations of subs() and map() use it, so you have
to implement either let_op(), or subs() and map().
You can, of course, also add your own new member functions. Remember
that the RTTI may be used to get information about what kinds of objects
you are dealing with (the position in the class hierarchy) and that you
can always extract the bare object from an ex by stripping the
ex off using the ex_to<mystring>(e) function when that
should become a need.
That's it. May the source be with you!
| [ << ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |