]> www.ginac.de Git - ginac.git/blobdiff - doc/tutorial/ginac.texi
degree(), ldegree(), coeff(), lcoeff(), tcoeff() and collect() can now
[ginac.git] / doc / tutorial / ginac.texi
index c85eadabf7611e3696b8e39b410ea15dc0e6a34c..d3a4d992beb56c90b726a9185f896a39874c8271 100644 (file)
@@ -23,7 +23,7 @@
 This is a tutorial that documents GiNaC @value{VERSION}, an open
 framework for symbolic computation within the C++ programming language.
 
-Copyright (C) 1999-2000 Johannes Gutenberg University Mainz, Germany
+Copyright (C) 1999-2001 Johannes Gutenberg University Mainz, Germany
 
 Permission is granted to make and distribute verbatim copies of
 this manual provided the copyright notice and this permission notice
@@ -52,7 +52,7 @@ notice identical to this one.
 
 @page
 @vskip 0pt plus 1filll
-Copyright @copyright{} 1999-2000 Johannes Gutenberg University Mainz, Germany
+Copyright @copyright{} 1999-2001 Johannes Gutenberg University Mainz, Germany
 @sp 2
 Permission is granted to make and distribute verbatim copies of
 this manual provided the copyright notice and this permission notice
@@ -135,7 +135,7 @@ the near future.
 
 @section License
 The GiNaC framework for symbolic computation within the C++ programming
-language is Copyright @copyright{} 1999-2000 Johannes Gutenberg
+language is Copyright @copyright{} 1999-2001 Johannes Gutenberg
 University Mainz, Germany.
 
 This program is free software; you can redistribute it and/or
@@ -180,6 +180,7 @@ pointless) bivariate polynomial with some large coefficients:
 
 @example
 #include <ginac/ginac.h>
+using namespace std;
 using namespace GiNaC;
 
 int main()
@@ -213,6 +214,7 @@ generates Hermite polynomials in a specified free variable.
 
 @example
 #include <ginac/ginac.h>
+using namespace std;
 using namespace GiNaC;
 
 ex HermitePoly(const symbol & x, int n)
@@ -661,6 +663,7 @@ meta-class for storing all mathematical objects.
 * Lists::                        Lists of expressions.
 * Mathematical functions::       Mathematical functions.
 * Relations::                    Equality, Inequality and all that.
+* Indexed objects::              Handling indexed quantities.
 @end menu
 
 
@@ -754,9 +757,10 @@ $\sqrt{2}$
 @item @code{lst} @tab Lists of expressions [@math{x}, @math{2*y}, @math{3+z}]
 @item @code{matrix} @tab @math{n}x@math{m} matrices of expressions
 @item @code{relational} @tab A relation like the identity @math{x}@code{==}@math{y}
-@item @code{color}, @code{coloridx} @tab Element and index of the @math{SU(3)} Lie-algebra
-@item @code{isospin} @tab Element of the @math{SU(2)} Lie-algebra
-@item @code{idx} @tab Index of a general tensor object
+@item @code{indexed} @tab Indexed object like @math{A_ij}
+@item @code{tensor} @tab Special tensor like the delta and metric tensors
+@item @code{idx} @tab Index of an indexed object
+@item @code{varidx} @tab Index with variance
 @end multitable
 @end cartouche
 
@@ -847,7 +851,7 @@ int main()
     // Trott's constant in scientific notation:
     numeric trott("1.0841015122311136151E-2");
     
-    cout << two*p << endl;  // floating point 6.283...
+    std::cout << two*p << std::endl;  // floating point 6.283...
 @}
 @end example
 
@@ -886,6 +890,7 @@ digits:
 
 @example
 #include <ginac/ginac.h>
+using namespace std;
 using namespace GiNaC;
 
 void foo()
@@ -936,6 +941,7 @@ some multiple of its denominator and test what comes out:
 
 @example
 #include <ginac/ginac.h>
+using namespace std;
 using namespace GiNaC;
 
 // some very important constants:
@@ -1040,21 +1046,16 @@ Simple polynomial expressions are written down in GiNaC pretty much like
 in other CAS or like expressions involving numerical variables in C.
 The necessary operators @code{+}, @code{-}, @code{*} and @code{/} have
 been overloaded to achieve this goal.  When you run the following
-program, the constructor for an object of type @code{mul} is
+code snippet, the constructor for an object of type @code{mul} is
 automatically called to hold the product of @code{a} and @code{b} and
 then the constructor for an object of type @code{add} is called to hold
 the sum of that @code{mul} object and the number one:
 
 @example
-#include <ginac/ginac.h>
-using namespace GiNaC;
-
-int main()
-@{
+    ...
     symbol a("a"), b("b");
     ex MyTerm = 1+a*b;
-    // ...
-@}
+    ...
 @end example
 
 @cindex @code{pow()}
@@ -1062,7 +1063,7 @@ For exponentiation, you have already seen the somewhat clumsy (though C-ish)
 statement @code{pow(x,2);} to represent @code{x} squared.  This direct
 construction is necessary since we cannot safely overload the constructor
 @code{^} in C++ to construct a @code{power} object.  If we did, it would
-have several counterintuitive effects:
+have several counterintuitive and undesired effects:
 
 @itemize @bullet
 @item
@@ -1171,37 +1172,27 @@ There are quite a number of useful functions hard-wired into GiNaC.  For
 instance, all trigonometric and hyperbolic functions are implemented
 (@xref{Built-in Functions}, for a complete list).
 
-These functions are all objects of class @code{function}.  They accept one
-or more expressions as arguments and return one expression.  If the arguments
-are not numerical, the evaluation of the function may be halted, as it
-does in the next example:
+These functions are all objects of class @code{function}.  They accept
+one or more expressions as arguments and return one expression.  If the
+arguments are not numerical, the evaluation of the function may be
+halted, as it does in the next example, showing how a function returns
+itself twice and finally an expression that may be really useful:
 
 @cindex Gamma function
 @cindex @code{subs()}
 @example
-#include <ginac/ginac.h>
-using namespace GiNaC;
-
-int main()
-@{
-    symbol x("x"), y("y");
-    
+    ...
+    symbol x("x"), y("y");    
     ex foo = x+y/2;
-    cout << "tgamma(" << foo << ") -> " << tgamma(foo) << endl;
+    cout << tgamma(foo) << endl;
+     // -> tgamma(x+(1/2)*y)
     ex bar = foo.subs(y==1);
-    cout << "tgamma(" << bar << ") -> " << tgamma(bar) << endl;
+    cout << tgamma(bar) << endl;
+     // -> tgamma(x+1/2)
     ex foobar = bar.subs(x==7);
-    cout << "tgamma(" << foobar << ") -> " << tgamma(foobar) << endl;
-@}
-@end example
-
-This program shows how the function returns itself twice and finally an
-expression that may be really useful:
-
-@example
-tgamma(x+(1/2)*y) -> tgamma(x+(1/2)*y)
-tgamma(x+1/2) -> tgamma(x+1/2)
-tgamma(15/2) -> (135135/128)*Pi^(1/2)
+    cout << tgamma(foobar) << endl;
+     // -> (135135/128)*Pi^(1/2)
+    ...
 @end example
 
 Besides evaluation most of these functions allow differentiation, series
@@ -1209,7 +1200,7 @@ expansion and so on.  Read the next chapter in order to learn more about
 this.
 
 
-@node Relations, Methods and Functions, Mathematical functions, Basic Concepts
+@node Relations, Indexed objects, Mathematical functions, Basic Concepts
 @c    node-name, next, previous, up
 @section Relations
 @cindex @code{relational} (class)
@@ -1236,7 +1227,537 @@ however, that @code{==} here does not perform any simplifications, hence
 @code{expand()} must be called explicitly.
 
 
-@node Methods and Functions, Information About Expressions, Relations, Top
+@node Indexed objects, Methods and Functions, Relations, Basic Concepts
+@c    node-name, next, previous, up
+@section Indexed objects
+
+GiNaC allows you to handle expressions containing general indexed objects in
+arbitrary spaces. It is also able to canonicalize and simplify such
+expressions and perform symbolic dummy index summations. There are a number
+of predefined indexed objects provided, like delta and metric tensors.
+
+There are few restrictions placed on indexed objects and their indices and
+it is easy to construct nonsense expressions, but our intention is to
+provide a general framework that allows you to implement algorithms with
+indexed quantities, getting in the way as little as possible.
+
+@cindex @code{idx} (class)
+@cindex @code{indexed} (class)
+@subsection Indexed quantities and their indices
+
+Indexed expressions in GiNaC are constructed of two special types of objects,
+@dfn{index objects} and @dfn{indexed objects}.
+
+@itemize @bullet
+
+@cindex contravariant
+@cindex covariant
+@cindex variance
+@item Index objects are of class @code{idx} or a subclass. Every index has
+a @dfn{value} and a @dfn{dimension} (which is the dimension of the space
+the index lives in) which can both be arbitrary expressions but are usually
+a number or a simple symbol. In addition, indices of class @code{varidx} have
+a @dfn{variance} (they can be co- or contravariant).
+
+@item Indexed objects are of class @code{indexed} or a subclass. They
+contain a @dfn{base expression} (which is the expression being indexed), and
+one or more indices.
+
+@end itemize
+
+@strong{Note:} when printing expressions, covariant indices and indices
+without variance are denoted @samp{.i} while contravariant indices are denoted
+@samp{~i}. In the following, we are going to use that notation in the text
+so instead of @math{A^i_jk} we will write @samp{A~i.j.k}. Index dimensions
+are not visible in the output.
+
+A simple example shall illustrate the concepts:
+
+@example
+#include <ginac/ginac.h>
+using namespace std;
+using namespace GiNaC;
+
+int main()
+@{
+    symbol i_sym("i"), j_sym("j");
+    idx i(i_sym, 3), j(j_sym, 3);
+
+    symbol A("A");
+    cout << indexed(A, i, j) << endl;
+     // -> A.i.j
+    ...
+@end example
+
+The @code{idx} constructor takes two arguments, the index value and the
+index dimension. First we define two index objects, @code{i} and @code{j},
+both with the numeric dimension 3. The value of the index @code{i} is the
+symbol @code{i_sym} (which prints as @samp{i}) and the value of the index
+@code{j} is the symbol @code{j_sym} (which prints as @samp{j}). Next we
+construct an expression containing one indexed object, @samp{A.i.j}. It has
+the symbol @code{A} as its base expression and the two indices @code{i} and
+@code{j}.
+
+Note the difference between the indices @code{i} and @code{j} which are of
+class @code{idx}, and the index values which are the sybols @code{i_sym}
+and @code{j_sym}. The indices of indexed objects cannot directly be symbols
+or numbers but must be index objects. For example, the following is not
+correct and will raise an exception:
+
+@example
+symbol i("i"), j("j");
+e = indexed(A, i, j); // ERROR: indices must be of type idx
+@end example
+
+You can have multiple indexed objects in an expression, index values can
+be numeric, and index dimensions symbolic:
+
+@example
+    ...
+    symbol B("B"), dim("dim");
+    cout << 4 * indexed(A, i)
+          + indexed(B, idx(j_sym, 4), idx(2, 3), idx(i_sym, dim)) << endl;
+     // -> B.j.2.i+4*A.i
+    ...
+@end example
+
+@code{B} has a 4-dimensional symbolic index @samp{k}, a 3-dimensional numeric
+index of value 2, and a symbolic index @samp{i} with the symbolic dimension
+@samp{dim}. Note that GiNaC doesn't automatically notify you that the free
+indices of @samp{A} and @samp{B} in the sum don't match (you have to call
+@code{simplify_indexed()} for that, see below).
+
+In fact, base expressions, index values and index dimensions can be
+arbitrary expressions:
+
+@example
+    ...
+    cout << indexed(A+B, idx(2*i_sym+1, dim/2)) << endl;
+     // -> (B+A).(1+2*i)
+    ...
+@end example
+
+It's also possible to construct nonsense like @samp{Pi.sin(x)}. You will not
+get an error message from this but you will probably not be able to do
+anything useful with it.
+
+@cindex @code{get_value()}
+@cindex @code{get_dimension()}
+The methods
+
+@example
+ex idx::get_value(void);
+ex idx::get_dimension(void);
+@end example
+
+return the value and dimension of an @code{idx} object. If you have an index
+in an expression, such as returned by calling @code{.op()} on an indexed
+object, you can get a reference to the @code{idx} object with the function
+@code{ex_to_idx()} on the expression.
+
+There are also the methods
+
+@example
+bool idx::is_numeric(void);
+bool idx::is_symbolic(void);
+bool idx::is_dim_numeric(void);
+bool idx::is_dim_symbolic(void);
+@end example
+
+for checking whether the value and dimension are numeric or symbolic
+(non-numeric). Using the @code{info()} method of an index (see @ref{Information
+About Expressions}) returns information about the index value.
+
+@cindex @code{varidx} (class)
+If you need co- and contravariant indices, use the @code{varidx} class:
+
+@example
+    ...
+    symbol mu_sym("mu"), nu_sym("nu");
+    varidx mu(mu_sym, 4), nu(nu_sym, 4); // default is contravariant ~mu, ~nu
+    varidx mu_co(mu_sym, 4, true);       // covariant index .mu
+
+    cout << indexed(A, mu, nu) << endl;
+     // -> A~mu~nu
+    cout << indexed(A, mu_co, nu) << endl;
+     // -> A.mu~nu
+    cout << indexed(A, mu.toggle_variance(), nu) << endl;
+     // -> A.mu~nu
+    ...
+@end example
+
+A @code{varidx} is an @code{idx} with an additional flag that marks it as
+co- or contravariant. The default is a contravariant (upper) index, but
+this can be overridden by supplying a third argument to the @code{varidx}
+constructor. The two methods
+
+@example
+bool varidx::is_covariant(void);
+bool varidx::is_contravariant(void);
+@end example
+
+allow you to check the variance of a @code{varidx} object (use @code{ex_to_varidx()}
+to get the object reference from an expression). There's also the very useful
+method
+
+@example
+ex varidx::toggle_variance(void);
+@end example
+
+which makes a new index with the same value and dimension but the opposite
+variance. By using it you only have to define the index once.
+
+@subsection Substituting indices
+
+@cindex @code{subs()}
+Sometimes you will want to substitute one symbolic index with another
+symbolic or numeric index, for example when calculating one specific element
+of a tensor expression. This is done with the @code{.subs()} method, as it
+is done for symbols (see @ref{Substituting Symbols}).
+
+You have two possibilities here. You can either substitute the whole index
+by another index or expression:
+
+@example
+    ...
+    ex e = indexed(A, mu_co);
+    cout << e << " becomes " << e.subs(mu_co == nu) << endl;
+     // -> A.mu becomes A~nu
+    cout << e << " becomes " << e.subs(mu_co == varidx(0, 4)) << endl;
+     // -> A.mu becomes A~0
+    cout << e << " becomes " << e.subs(mu_co == 0) << endl;
+     // -> A.mu becomes A.0
+    ...
+@end example
+
+The third example shows that trying to replace an index with something that
+is not an index will substitute the index value instead.
+
+Alternatively, you can substitute the @emph{symbol} of a symbolic index by
+another expression:
+
+@example
+    ...
+    ex e = indexed(A, mu_co);
+    cout << e << " becomes " << e.subs(mu_sym == nu_sym) << endl;
+     // -> A.mu becomes A.nu
+    cout << e << " becomes " << e.subs(mu_sym == 0) << endl;
+     // -> A.mu becomes A.0
+    ...
+@end example
+
+As you see, with the second method only the value of the index will get
+substituted. Its other properties, including its dimension, remain unchanged.
+If you want to change the dimension of an index you have to substitute the
+whole index by another one with the new dimension.
+
+Finally, substituting the base expression of an indexed object works as
+expected:
+
+@example
+    ...
+    ex e = indexed(A, mu_co);
+    cout << e << " becomes " << e.subs(A == A+B) << endl;
+     // -> A.mu becomes (B+A).mu
+    ...
+@end example
+
+@subsection Symmetries
+
+Indexed objects can be declared as being totally symmetric or antisymmetric
+with respect to their indices. In this case, GiNaC will automatically bring
+the indices into a canonical order which allows for some immediate
+simplifications:
+
+@example
+    ...
+    cout << indexed(A, indexed::symmetric, i, j)
+          + indexed(A, indexed::symmetric, j, i) << endl;
+     // -> 2*A.j.i
+    cout << indexed(B, indexed::antisymmetric, i, j)
+          + indexed(B, indexed::antisymmetric, j, j) << endl;
+     // -> -B.j.i
+    cout << indexed(B, indexed::antisymmetric, i, j)
+          + indexed(B, indexed::antisymmetric, j, i) << endl;
+     // -> 0
+    ...
+@end example
+
+@cindex @code{get_free_indices()}
+@cindex Dummy index
+@subsection Dummy indices
+
+GiNaC treats certain symbolic index pairs as @dfn{dummy indices} meaning
+that a summation over the index range is implied. Symbolic indices which are
+not dummy indices are called @dfn{free indices}. Numeric indices are neither
+dummy nor free indices.
+
+To be recognized as a dummy index pair, the two indices must be of the same
+class and dimension and their value must be the same single symbol (an index
+like @samp{2*n+1} is never a dummy index). If the indices are of class
+@code{varidx}, they must also be of opposite variance.
+
+The method @code{.get_free_indices()} returns a vector containing the free
+indices of an expression. It also checks that the free indices of the terms
+of a sum are consistent:
+
+@example
+@{
+    symbol A("A"), B("B"), C("C");
+
+    symbol i_sym("i"), j_sym("j"), k_sym("k"), l_sym("l");
+    idx i(i_sym, 3), j(j_sym, 3), k(k_sym, 3), l(l_sym, 3);
+
+    ex e = indexed(A, i, j) * indexed(B, j, k) + indexed(C, k, l, i, l);
+    cout << exprseq(e.get_free_indices()) << endl;
+     // -> (.i,.k)
+     // 'j' and 'l' are dummy indices
+
+    symbol mu_sym("mu"), nu_sym("nu"), rho_sym("rho"), sigma_sym("sigma");
+    varidx mu(mu_sym, 4), nu(nu_sym, 4), rho(rho_sym, 4), sigma(sigma_sym, 4);
+
+    e = indexed(A, mu, nu) * indexed(B, nu.toggle_variance(), rho)
+      + indexed(C, mu, sigma, rho, sigma.toggle_variance());
+    cout << exprseq(e.get_free_indices()) << endl;
+     // -> (~mu,~rho)
+     // 'nu' is a dummy index, but 'sigma' is not
+
+    e = indexed(A, mu, mu);
+    cout << exprseq(e.get_free_indices()) << endl;
+     // -> (~mu)
+     // 'mu' is not a dummy index because it appears twice with the same
+     // variance
+
+    e = indexed(A, mu, nu) + 42;
+    cout << exprseq(e.get_free_indices()) << endl; // ERROR
+     // this will throw an exception:
+     // "add::get_free_indices: inconsistent indices in sum"
+@}
+@end example
+
+@cindex @code{simplify_indexed()}
+@subsection Simplifying indexed expressions
+
+In addition to the few automatic simplifications that GiNaC performs on
+indexed expressions (such as re-ordering the indices of symmetric tensors
+and calculating traces and convolutions of matrices and predefined tensors)
+there is the method
+
+@example
+ex ex::simplify_indexed(void);
+ex ex::simplify_indexed(const scalar_products & sp);
+@end example
+
+that performs some more expensive operations:
+
+@itemize
+@item it checks the consistency of free indices in sums in the same way
+  @code{get_free_indices()} does
+@item it (symbolically) calculates all possible dummy index summations/contractions
+  with the predefined tensors (this will be explained in more detail in the
+  next section)
+@item as a special case of dummy index summation, it can replace scalar products
+  of two tensors with a user-defined value
+@end itemize
+
+The last point is done with the help of the @code{scalar_products} class
+which is used to store scalar products with known values (this is not an
+arithmetic class, you just pass it to @code{simplify_indexed()}):
+
+@example
+@{
+    symbol A("A"), B("B"), C("C"), i_sym("i");
+    idx i(i_sym, 3);
+
+    scalar_products sp;
+    sp.add(A, B, 0); // A and B are orthogonal
+    sp.add(A, C, 0); // A and C are orthogonal
+    sp.add(A, A, 4); // A^2 = 4 (A has length 2)
+
+    e = indexed(A + B, i) * indexed(A + C, i);
+    cout << e << endl;
+     // -> (B+A).i*(A+C).i
+
+    cout << e.expand(expand_options::expand_indexed).simplify_indexed(sp)
+         << endl;
+     // -> 4+C.i*B.i
+@}
+@end example
+
+The @code{scalar_products} object @code{sp} acts as a storage for the
+scalar products added to it with the @code{.add()} method. This method
+takes three arguments: the two expressions of which the scalar product is
+taken, and the expression to replace it with. After @code{sp.add(A, B, 0)},
+@code{simplify_indexed()} will replace all scalar products of indexed
+objects that have the symbols @code{A} and @code{B} as base expressions
+with the single value 0. The number, type and dimension of the indices
+doesn't matter; @samp{A~mu~nu*B.mu.nu} would also be replaced by 0.
+
+@cindex @code{expand()}
+The example above also illustrates a feature of the @code{expand()} method:
+if passed the @code{expand_indexed} option it will distribute indices
+over sums, so @samp{(A+B).i} becomes @samp{A.i+B.i}.
+
+@cindex @code{tensor} (class)
+@subsection Predefined tensors
+
+Some frequently used special tensors such as the delta, epsilon and metric
+tensors are predefined in GiNaC. They have special properties when
+contracted with other tensor expressions and some of them have constant
+matrix representations (they will evaluate to a number when numeric
+indices are specified).
+
+@cindex @code{delta_tensor()}
+@subsubsection Delta tensor
+
+The delta tensor takes two indices, is symmetric and has the matrix
+representation @code{diag(1,1,1,...)}. It is constructed by the function
+@code{delta_tensor()}:
+
+@example
+@{
+    symbol A("A"), B("B");
+
+    idx i(symbol("i"), 3), j(symbol("j"), 3),
+        k(symbol("k"), 3), l(symbol("l"), 3);
+
+    ex e = indexed(A, i, j) * indexed(B, k, l)
+         * delta_tensor(i, k) * delta_tensor(j, l) << endl;
+    cout << e.simplify_indexed() << endl;
+     // -> B.i.j*A.i.j
+
+    cout << delta_tensor(i, i) << endl;
+     // -> 3
+@}
+@end example
+
+@cindex @code{metric_tensor()}
+@subsubsection General metric tensor
+
+The function @code{metric_tensor()} creates a general symmetric metric
+tensor with two indices that can be used to raise/lower tensor indices. The
+metric tensor is denoted as @samp{g} in the output and if its indices are of
+mixed variance it is automatically replaced by a delta tensor:
+
+@example
+@{
+    symbol A("A");
+
+    varidx mu(symbol("mu"), 4), nu(symbol("nu"), 4), rho(symbol("rho"), 4);
+
+    ex e = metric_tensor(mu, nu) * indexed(A, nu.toggle_variance(), rho);
+    cout << e.simplify_indexed() << endl;
+     // -> A~mu~rho
+
+    e = delta_tensor(mu, nu.toggle_variance()) * metric_tensor(nu, rho);
+    cout << e.simplify_indexed() << endl;
+     // -> g~mu~rho
+
+    e = metric_tensor(mu.toggle_variance(), nu.toggle_variance())
+      * metric_tensor(nu, rho);
+    cout << e.simplify_indexed() << endl;
+     // -> delta.mu~rho
+
+    e = metric_tensor(nu.toggle_variance(), rho.toggle_variance())
+      * metric_tensor(mu, nu) * (delta_tensor(mu.toggle_variance(), rho)
+        + indexed(A, mu.toggle_variance(), rho));
+    cout << e.simplify_indexed() << endl;
+     // -> 4+A.rho~rho
+@}
+@end example
+
+@cindex @code{lorentz_g()}
+@subsubsection Minkowski metric tensor
+
+The Minkowski metric tensor is a special metric tensor with a constant
+matrix representation which is either @code{diag(1, -1, -1, ...)} (negative
+signature, the default) or @code{diag(-1, 1, 1, ...)} (positive signature).
+It is created with the function @code{lorentz_g()} (although it is output as
+@samp{eta}):
+
+@example
+@{
+    varidx mu(symbol("mu"), 4);
+
+    e = delta_tensor(varidx(0, 4), mu.toggle_variance())
+      * lorentz_g(mu, varidx(0, 4));       // negative signature
+    cout << e.simplify_indexed() << endl;
+     // -> 1
+
+    e = delta_tensor(varidx(0, 4), mu.toggle_variance())
+      * lorentz_g(mu, varidx(0, 4), true); // positive signature
+    cout << e.simplify_indexed() << endl;
+     // -> -1
+@}
+@end example
+
+@subsubsection Epsilon tensor
+
+The epsilon tensor is totally antisymmetric, its number of indices is equal
+to the dimension of the index space (the indices must all be of the same
+numeric dimension), and @samp{eps.1.2.3...} (resp. @samp{eps~0~1~2...}) is
+defined to be 1. Its behaviour with indices that have a variance also
+depends on the signature of the metric. Epsilon tensors are output as
+@samp{eps}.
+
+There are three functions defined to create epsilon tensors in 2, 3 and 4
+dimensions:
+
+@example
+ex epsilon_tensor(const ex & i1, const ex & i2);
+ex epsilon_tensor(const ex & i1, const ex & i2, const ex & i3);
+ex lorentz_eps(const ex & i1, const ex & i2, const ex & i3, const ex & i4, bool pos_sig = false);
+@end example
+
+The first two functions create an epsilon tensor in 2 or 3 Euclidean
+dimensions, the last function creates an epsilon tensor in a 4-dimensional
+Minkowski space (the last @code{bool} argument specifies whether the metric
+has negative or positive signature, as in the case of the Minkowski metric
+tensor).
+
+@subsection Linear algebra
+
+The @code{matrix} class can be used with indices to do some simple linear
+algebra (linear combinations and products of vectors and matrices, traces
+and scalar products):
+
+@example
+@{
+    idx i(symbol("i"), 2), j(symbol("j"), 2);
+    symbol x("x"), y("y");
+
+    matrix A(2, 2, lst(1, 2, 3, 4)), X(2, 1, lst(x, y));
+
+    cout << indexed(A, i, i) << endl;
+     // -> 5
+
+    ex e = indexed(A, i, j) * indexed(X, j);
+    cout << e.simplify_indexed() << endl;
+     // -> [[ [[2*y+x]], [[4*y+3*x]] ]].i
+
+    e = indexed(A, i, j) * indexed(X, i) + indexed(X, j) * 2;
+    cout << e.simplify_indexed() << endl;
+     // -> [[ [[3*y+3*x,6*y+2*x]] ]].j
+@}
+@end example
+
+You can of course obtain the same results with the @code{matrix::add()},
+@code{matrix::mul()} and @code{matrix::trace()} methods but with indices you
+don't have to worry about transposing matrices.
+
+Matrix indices always start at 0 and their dimension must match the number
+of rows/columns of the matrix. Matrices with one row or one column are
+vectors and can have one or two indices (it doesn't matter whether it's a
+row or a column vector). Other matrices must have two indices.
+
+You should be careful when using indices with variance on matrices. GiNaC
+doesn't look at the variance and doesn't know that @samp{F~mu~nu} and
+@samp{F.mu.nu} are different matrices. In this case you should use only
+one form for @samp{F} and explicitly multiply it with a matrix representation
+of the metric tensor.
+
+
+@node Methods and Functions, Information About Expressions, Indexed objects, Top
 @c    node-name, next, previous, up
 @chapter Methods and Functions
 @cindex polynomial
@@ -1249,16 +1770,10 @@ alternatively call it in a functional way as shown in the simple
 example:
 
 @example
-#include <ginac/ginac.h>
-using namespace GiNaC;
-
-int main()
-@{
-    ex x = numeric(1.0);
-    
-    cout << "As method:   " << sin(x).evalf() << endl;
-    cout << "As function: " << evalf(sin(x)) << endl;
-@}
+    ...
+    cout << "As method:   " << sin(1).evalf() << endl;
+    cout << "As function: " << evalf(sin(1)) << endl;
+    ...
 @end example
 
 @cindex @code{subs()}
@@ -1300,6 +1815,9 @@ avoided.
 
 @subsection Checking expression types
 @cindex @code{is_ex_of_type()}
+@cindex @code{ex_to_numeric()}
+@cindex @code{ex_to_@dots{}}
+@cindex @code{Converting ex to other classes}
 @cindex @code{info()}
 
 Sometimes it's useful to check whether a given expression is a plain number,
@@ -1311,6 +1829,20 @@ bool is_ex_of_type(const ex & e, TYPENAME t);
 bool ex::info(unsigned flag);
 @end example
 
+When the test made by @code{is_ex_of_type()} returns true, it is safe to
+call one of the functions @code{ex_to_@dots{}}, where @code{@dots{}} is
+one of the class names (@xref{The Class Hierarchy}, for a list of all
+classes). For example, assuming @code{e} is an @code{ex}:
+
+@example
+@{
+    @dots{}
+    if (is_ex_of_type(e, numeric))
+        numeric n = ex_to_numeric(e);
+    @dots{}
+@}
+@end example
+
 @code{is_ex_of_type()} allows you to check whether the top-level object of
 an expression @samp{e} is an instance of the GiNaC class @samp{t}
 (@xref{The Class Hierarchy}, for a list of all classes). This is most useful,
@@ -1397,7 +1929,9 @@ table:
 @item @code{crational_polynomial}
 @tab @dots{}a polynomial with (possibly complex) rational coefficients (such as @math{2/3+7/2*I})
 @item @code{rational_function}
-@tab @dots{}a rational function
+@tab @dots{}a rational function (@math{x+y}, @math{z/(x+y)})
+@item @code{algebraic}
+@tab @dots{}an algebraic object (@math{sqrt(2)}, @math{sqrt(x)-1})
 @end multitable
 @end cartouche
 
@@ -1421,7 +1955,8 @@ for accessing the subexpressions in the container-like GiNaC classes like
 determines the number of subexpressions (@samp{operands}) contained, while
 @code{op()} returns the @code{i}-th (0..@code{nops()-1}) subexpression.
 In the case of a @code{power} object, @code{op(0)} will return the basis
-and @code{op(1)} the exponent.
+and @code{op(1)} the exponent. For @code{indexed} objects, @code{op(0)}
+is the base expression and @code{op(i)}, @math{i>0} are the indices.
 
 The left-hand and right-hand side expressions of objects of class
 @code{relational} (and only of these) can also be accessed with the methods
@@ -1548,7 +2083,7 @@ being polynomials in the remaining variables.  The method
 @code{collect()} accomplishes this task:
 
 @example
-ex ex::collect(const symbol & s);
+ex ex::collect(const ex & s);
 @end example
 
 Note that the original polynomial needs to be in expanded form in order
@@ -1563,8 +2098,8 @@ The degree and low degree of a polynomial can be obtained using the two
 methods
 
 @example
-int ex::degree(const symbol & s);
-int ex::ldegree(const symbol & s);
+int ex::degree(const ex & s);
+int ex::ldegree(const ex & s);
 @end example
 
 which also work reliably on non-expanded input polynomials (they even work
@@ -1572,14 +2107,14 @@ on rational functions, returning the asymptotic degree). To extract
 a coefficient with a certain power from an expanded polynomial you use
 
 @example
-ex ex::coeff(const symbol & s, int n);
+ex ex::coeff(const ex & s, int n);
 @end example
 
 You can also obtain the leading and trailing coefficients with the methods
 
 @example
-ex ex::lcoeff(const symbol & s);
-ex ex::tcoeff(const symbol & s);
+ex ex::lcoeff(const ex & s);
+ex ex::tcoeff(const ex & s);
 @end example
 
 which are equivalent to @code{coeff(s, degree(s))} and @code{coeff(s, ldegree(s))},
@@ -1590,6 +2125,7 @@ polynomial is analyzed:
 
 @example
 #include <ginac/ginac.h>
+using namespace std;
 using namespace GiNaC;
 
 int main()
@@ -1722,6 +2258,39 @@ int main()
 @end example
 
 
+@subsection Square-free decomposition
+@cindex square-free decomposition
+@cindex factorization
+@cindex @code{sqrfree()}
+
+GiNaC still lacks proper factorization support.  Some form of
+factorization is, however, easily implemented by noting that factors
+appearing in a polynomial with power two or more also appear in the
+derivative and hence can easily be found by computing the GCD of the
+original polynomial and its derivatives.  Any system has an interface
+for this so called square-free factorization.  So we provide one, too:
+@example
+ex sqrfree(const ex & a, const lst & l = lst());
+@end example
+Here is an example that by the way illustrates how the result may depend
+on the order of differentiation:
+@example
+    ...
+    symbol x("x"), y("y");
+    ex BiVarPol = expand(pow(x-2*y*x,3) * pow(x+y,2) * (x-y));
+
+    cout << sqrfree(BiVarPol, lst(x,y)) << endl;
+     // -> (y+x)^2*(-1+6*y+8*y^3-12*y^2)*(y-x)*x^3
+
+    cout << sqrfree(BiVarPol, lst(y,x)) << endl;
+     // -> (1-2*y)^3*(y+x)^2*(-y+x)*x^3
+
+    cout << sqrfree(BiVarPol) << endl;
+     // -> depending on luck, any of the above
+    ...
+@end example
+
+
 @node Rational Expressions, Symbolic Differentiation, Polynomial Arithmetic, Methods and Functions
 @c    node-name, next, previous, up
 @section Rational expressions
@@ -1731,7 +2300,7 @@ int main()
 @cindex simplification
 @cindex temporary replacement
 
-Some basic from of simplification of expressions is called for frequently.
+Some basic form of simplification of expressions is called for frequently.
 GiNaC provides the method @code{.normal()}, which converts a rational function
 into an equivalent rational function of the form @samp{numerator/denominator}
 where numerator and denominator are coprime.  If the input expression is already
@@ -1757,8 +2326,8 @@ int main()
     symbol x("x");
     ex t1 = (pow(x,2) + 2*x + 1)/(x + 1);
     ex t2 = (pow(sin(x),2) + 2*sin(x) + 1)/(sin(x) + 1);
-    cout << "t1 is " << t1.normal() << endl;
-    cout << "t2 is " << t2.normal() << endl;
+    std::cout << "t1 is " << t1.normal() << std::endl;
+    std::cout << "t2 is " << t2.normal() << std::endl;
 @}
 @end example
 
@@ -1877,7 +2446,7 @@ ex EulerNumber(unsigned n)
 int main()
 @{
     for (unsigned i=0; i<11; i+=2)
-        cout << EulerNumber(i) << endl;
+        std::cout << EulerNumber(i) << std::endl;
     return 0;
 @}
 @end example
@@ -1907,6 +2476,7 @@ term).  A sample application from special relativity could read:
 
 @example
 #include <ginac/ginac.h>
+using namespace std;
 using namespace GiNaC;
 
 int main()
@@ -1964,6 +2534,8 @@ ex mechain_pi(int degr)
 
 int main()
 @{
+    using std::cout;  // just for fun, another way of...
+    using std::endl;  // ...dealing with this namespace std.
     ex pi_frac;
     for (int i=2; i<12; i+=2) @{
         pi_frac = mechain_pi(i);
@@ -2039,6 +2611,8 @@ GiNaC contains the following predefined mathematical functions:
 @tab exponential function
 @item @code{log(x)}
 @tab natural logarithm
+@item @code{Li2(x)}
+@tab Dilogarithm
 @item @code{zeta(x)}
 @tab Riemann's zeta function
 @item @code{zeta(n, x)}
@@ -2066,10 +2640,19 @@ GiNaC contains the following predefined mathematical functions:
 
 @cindex branch cut
 For functions that have a branch cut in the complex plane GiNaC follows
-the conventions for C++ as defined in the ANSI standard.  In particular:
-the natural logarithm (@code{log}) and the square root (@code{sqrt})
-both have their branch cuts running along the negative real axis where
-the points on the axis itself belong to the upper part.
+the conventions for C++ as defined in the ANSI standard as far as
+possible.  In particular: the natural logarithm (@code{log}) and the
+square root (@code{sqrt}) both have their branch cuts running along the
+negative real axis where the points on the axis itself belong to the
+upper part (i.e. continuous with quadrant II).  The inverse
+trigonometric and hyperbolic functions are not defined for complex
+arguments by the C++ standard, however.  In GiNaC we follow the
+conventions used by CLN, which in turn follow the carefully designed
+definitions in the Common Lisp standard.  It should be noted that this
+convention is identical to the one used by the C99 standard and by most
+serious CAS.  It is to be expected that future revisions of the C++
+standard incorporate these functions in the complex domain in a manner
+compatible with C99.
 
 
 @node Input/Output, Extending GiNaC, Built-in Functions, Methods and Functions
@@ -2185,6 +2768,7 @@ With this constructor, it's also easy to implement interactive GiNaC programs:
 #include <string>
 #include <stdexcept>
 #include <ginac/ginac.h>
+using namespace std;
 using namespace GiNaC;
 
 int main()
@@ -2216,8 +2800,9 @@ of class @code{archive} and archive expressions in it, giving each
 expression a unique name:
 
 @example
-#include <ginac/ginac.h>
 #include <fstream>
+using namespace std;
+#include <ginac/ginac.h>
 using namespace GiNaC;
 
 int main()
@@ -2307,6 +2892,7 @@ authors---they will happily incorporate them into future versions.
 @menu
 * What does not belong into GiNaC::  What to avoid.
 * Symbolic functions::               Implementing symbolic functions.
+* Adding classes::                   Defining new algebraic classes.
 @end menu
 
 
@@ -2334,7 +2920,7 @@ inefficient.  For this purpose, the underlying foundation classes
 provided by @acronym{CLN} are much better suited.
 
 
-@node Symbolic functions, A Comparison With Other CAS, What does not belong into GiNaC, Extending GiNaC
+@node Symbolic functions, Adding classes, What does not belong into GiNaC, Extending GiNaC
 @c    node-name, next, previous, up
 @section Symbolic functions
 
@@ -2400,7 +2986,7 @@ enough to know how to differentiate.  But if the function you want to
 implement does have a pole somewhere in the complex plane, you need to
 write another method for Laurent expansion around that point.
 
-Now that all the ingrediences for @code{cos} have been set up, we need
+Now that all the ingredients for @code{cos} have been set up, we need
 to tell the system about it.  This is done by a macro and we are not
 going to descibe how it expands, please consult your preprocessor if you
 are curious:
@@ -2438,10 +3024,502 @@ mechanisms.  Please, have a look at the real implementation in GiNaC.
 assure you that functions are GiNaC's most macro-intense classes.  We
 have done our best to avoid macros where we can.)
 
+
+@node Adding classes, A Comparison With Other CAS, Symbolic functions, Extending GiNaC
+@c    node-name, next, previous, up
+@section Adding classes
+
+If you are doing some very specialized things with GiNaC you may find that
+you have to implement your own algebraic classes to fit your needs. This
+section will explain how to do this 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.
+
+@subsection GiNaC's run-time type information system
+
+@cindex hierarchy of classes
+@cindex RTTI
+All algebraic classes (that is, all classes that can appear in expressions)
+in GiNaC are direct or indirect subclasses of the class @code{basic}. So a
+@code{basic *} (which is essentially what an @code{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 @code{basic *} really is.
+Also, for the unarchiving of expressions it must be possible to find the
+@code{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 @file{<typeinfo>}) but for efficiency reasons GiNaC
+implements its own, simpler RTTI.
+
+The RTTI in GiNaC is based on two mechanisms:
+
+@itemize @bullet
+
+@item
+The @code{basic} class declares a member variable @code{tinfo_key} which
+holds an unsigned integer that identifies the object's class. These numbers
+are defined in the @file{tinfos.h} header file for the built-in GiNaC
+classes. They all start with @code{TINFO_}.
+
+@item
+By means of some clever tricks with static members, GiNaC maintains a list
+of information for all classes derived from @code{basic}. The information
+available includes the class names, the @code{tinfo_key}s, and pointers
+to the unarchiving functions. This class registry is defined in the
+@file{registrar.h} header file.
+
+@end itemize
+
+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 automatic) but don't worry, most of the work is simplified by
+macros.
+
+@subsection A minimalistic example
+
+Now we will start implementing a new class @code{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
+@code{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:
+
+@example
+#include <iostream>
+#include <string>   
+#include <stdexcept>
+using namespace std;
+
+#include <ginac/ginac.h>
+using namespace GiNaC;
+@end example
+
+The first thing we have to do is to define a @code{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 @file{tinfos.h} are modeled somewhat after the class hierarchy
+which is not a requirement but we are going to stick with this scheme:
+
+@example
+const unsigned TINFO_mystring = 0x42420001U;
+@end example
+
+Now we can write down the class declaration. The class stores a C++
+@code{string} and the user shall be able to construct a @code{mystring}
+object from a C or C++ string:
+
+@example
+class mystring : public basic
+@{
+    GINAC_DECLARE_REGISTERED_CLASS(mystring, basic)
+  
+public:
+    mystring(const string &s);
+    mystring(const char *s);
+
+private:
+    string str;
+@};
+
+GIANC_IMPLEMENT_REGISTERED_CLASS(mystring, basic)
+@end example
+
+The @code{GINAC_DECLARE_REGISTERED_CLASS} and @code{GINAC_IMPLEMENT_REGISTERED_CLASS}
+macros are defined in @file{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 @code{GINAC_DECLARE_REGISTERED_CLASS} should be
+the first line after the opening brace of the class definition. The
+@code{GINAC_IMPLEMENT_REGISTERED_CLASS} may appear anywhere else in the
+source (at global scope, of course, not inside a function).
+
+@code{GINAC_DECLARE_REGISTERED_CLASS} contains, among other things the
+declarations of the default and copy constructor, the destructor, the
+assignment operator and a couple of other functions that are required. It
+also defines a type @code{inherited} which refers to the superclass so you
+don't have to modify your code every time you shuffle around the class
+hierarchy. @code{GINAC_IMPLEMENT_REGISTERED_CLASS} implements the copy
+constructor, the destructor and the assignment operator.
+
+Now there are nine member functions we have to implement to get a working
+class:
+
+@itemize
+
+@item
+@code{mystring()}, the default constructor.
+
+@item
+@code{void destroy(bool call_parent)}, which is used in the destructor and the
+assignment operator to free dynamically allocated members. The @code{call_parent}
+specifies whether the @code{destroy()} function of the superclass is to be
+called also.
+
+@item
+@code{void copy(const mystring &other)}, which is used in the copy constructor
+and assignment operator to copy the member variables over from another
+object of the same class.
+
+@item
+@code{void archive(archive_node &n)}, the archiving function. This stores all
+information needed to reconstruct an object of this class inside an
+@code{archive_node}.
+
+@item
+@code{mystring(const archive_node &n, const lst &sym_lst)}, the unarchiving
+constructor. This constructs an instance of the class from the information
+found in an @code{archive_node}.
+
+@item
+@code{ex unarchive(const archive_node &n, const lst &sym_lst)}, the static
+unarchiving function. It constructs a new instance by calling the unarchiving
+constructor.
+
+@item
+@code{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 @code{other}
+object. If it returns 0, the objects are considered equal.
+@strong{Note:} This has nothing to do with the (numeric) ordering
+relationship expressed by @code{<}, @code{>=} etc (which cannot be defined
+for non-numeric classes). For example, @code{numeric(1).compare_same_type(numeric(2))}
+may return +1 even though 1 is clearly smaller than 2. Every GiNaC class
+must provide a @code{compare_same_type()} function, even those representing
+objects for which no reasonable algebraic ordering relationship can be
+defined.
+
+@item
+And, of course, @code{mystring(const string &s)} and @code{mystring(const char *s)}
+which are the two constructors we declared.
+
+@end itemize
+
+Let's proceed step-by-step. The default constructor looks like this:
+
+@example
+mystring::mystring() : inherited(TINFO_mystring)
+@{
+    // dynamically allocate resources here if required
+@}
+@end example
+
+The golden rule is that in all constructors you have to set the
+@code{tinfo_key} member to the @code{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 @code{basic} class provides
+a constructor that takes a @code{tinfo_key} value, which we are using here
+(remember that in our case @code{inherited = basic}). If the superclass
+didn't have such a constructor, we would have to set the @code{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 @code{str}
+member gets set to an empty string automatically). The constructor(s) are of
+course also the right place to allocate any dynamic resources you require.
+
+Next, the @code{destroy()} function:
+
+@example
+void mystring::destroy(bool call_parent)
+@{
+    // free dynamically allocated resources here if required
+    if (call_parent)
+        inherited::destroy(call_parent);
+@}
+@end example
+
+This function is where we free all dynamically allocated resources. We don't
+have any so we're not doing anything here, but if we had, for example, used
+a C-style @code{char *} to store our string, this would be the place to
+@code{delete[]} the string storage. If @code{call_parent} is true, we have
+to call the @code{destroy()} function of the superclass after we're done
+(to mimic C++'s automatic invocation of superclass destructors where
+@code{destroy()} is called from outside a destructor).
+
+The @code{copy()} function just copies over the member variables from
+another object:
+
+@example
+void mystring::copy(const mystring &other)
+@{
+    inherited::copy(other);
+    str = other.str;
+@}
+@end example
+
+We can simply overwrite the member variables here. There's no need to worry
+about dynamically allocated storage. The assignment operator (which is
+automatically defined by @code{GINAC_IMPLEMENT_REGISTERED_CLASS}, as you
+recall) calls @code{destroy()} before it calls @code{copy()}. You have to
+explicitly call the @code{copy()} function of the superclass here so
+all the member variables will get copied.
+
+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:
+
+@example
+void mystring::archive(archive_node &n) const
+@{
+    inherited::archive(n);
+    n.add_string("string", str);
+@}
+@end example
+
+The only thing that is really required is calling the @code{archive()}
+function of the superclass. Optionally, you can store all information you
+deem necessary for representing the object into the passed
+@code{archive_node}. We are just storing our string here. For more
+information on how the archiving works, consult the @file{archive.h} header
+file.
+
+The unarchiving constructor is basically the inverse of the archiving
+function:
+
+@example
+mystring::mystring(const archive_node &n, const lst &sym_lst) : inherited(n, sym_lst)
+@{
+    n.find_string("string", str);
+@}
+@end example
+
+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 @code{tinfo_key} here because it is done automatically
+by the unarchiving constructor of the @code{basic} class.
+
+Finally, the unarchiving function:
+
+@example
+ex mystring::unarchive(const archive_node &n, const lst &sym_lst)
+@{
+    return (new mystring(n, sym_lst))->setflag(status_flags::dynallocated);
+@}
+@end example
+
+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 @code{archive_node}s to global
+objects) you don't need a different implementation. For those who are
+interested: setting the @code{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 @code{compare_same_type()} function uses a provided function to compare
+the string members:
+
+@example
+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;
+@}
+@end example
+
+Although this function takes a @code{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 @math{A-B=0}), so you should compare
+all relevant member variables.
+
+Now the only thing missing is our two new constructors:
+
+@example
+mystring::mystring(const string &s) : inherited(TINFO_mystring), str(s)
+@{
+    // dynamically allocate resources here if required
+@}
+
+mystring::mystring(const char *s) : inherited(TINFO_mystring), str(s)
+@{
+    // dynamically allocate resources here if required
+@}
+@end example
+
+No surprises here. We set the @code{str} member from the argument and
+remember to pass the right @code{tinfo_key} to the @code{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:
+
+@example
+ex e = mystring("Hello, world!");
+cout << is_ex_of_type(e, mystring) << endl;
+ // -> 1 (true)
+
+cout << e.bp->class_name() << endl;
+ // -> mystring
+@end example
+
+Obviously it does. Let's see what the expression @code{e} looks like:
+
+@example
+cout << e << endl;
+ // -> [mystring object]
+@end example
+
+Hm, not exactly what we expect, but of course the @code{mystring} class
+doesn't yet know how to print itself. This is done in the @code{print()}
+member function. Let's say that we wanted to print the string surrounded
+by double quotes:
+
+@example
+class mystring : public basic
+@{
+    ...
+public:
+    void print(ostream &os, unsigned upper_precedence) const;
+    ...
+@};
+
+void mystring::print(ostream &os, unsigned upper_precedence) const
+@{
+    os << '\"' << str << '\"';
+@}
+@end example
+
+The @code{upper_precedence} argument is only required for container classes
+to correctly parenthesize the output. Let's try again to print the expression:
+
+@example
+cout << e << endl;
+ // -> "Hello, world!"
+@end example
+
+Much better. The @code{mystring} class can be used in arbitrary expressions:
+
+@example
+e += mystring("GiNaC rulez"); 
+cout << e << endl;
+ // -> "GiNaC rulez"+"Hello, world!"
+@end example
+
+(note that GiNaC's automatic term reordering is in effect here), or even
+
+@example
+e = pow(mystring("One string"), 2*sin(Pi-mystring("Another string")));
+cout << e << endl;
+ // -> "One string"^(2*sin(-"Another string"+Pi))
+@end example
+
+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:
+
+@example
+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
+@end example
+
+There's no way to, for example, make GiNaC's @code{add} class perform string
+concatenation. You would have to implement this yourself.
+
+@subsection Automatic evaluation
+
+@cindex @code{hold()}
+@cindex evaluation
+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 @code{eval()}. Let's say that
+we wanted all strings automatically converted to lowercase with
+non-alphabetic characters stripped, and empty strings removed:
+
+@example
+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 _ex0();
+    else
+        return mystring(new_str).hold();
+@}
+@end example
+
+The @code{level} argument is used to limit the recursion depth of the
+evaluation. We don't have any subexpressions in the @code{mystring} class
+so we are not concerned with this. If we had, we would call the @code{eval()}
+functions of the subexpressions with @code{level - 1} as the argument if
+@code{level != 1}. The @code{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
+@code{return this->hold();}.
+
+Let's confirm that it works:
+
+@example
+ex e = mystring("Hello, world!") + mystring("!?#");
+cout << e << endl;
+ // -> "helloworld"
+
+e = mystring("Wow!") + mystring("WOW") + mystring(" W ** o ** W");  
+cout << e << endl;
+ // -> 3*"wow"
+@end example
+
+@subsection Other member functions
+
+We have implemented only a small set of member functions to make the class
+work in the GiNaC framework. For a real algebraic class, there are probably
+some more functions that you will want to re-implement, such as
+@code{evalf()}, @code{series()} or @code{op()}. Have a look at @file{basic.h}
+or the header file of the class you want to make a subclass of to see
+what's there. You can, of course, also add your own new member functions.
+In this case you will probably want to define a little helper function like
+
+@example
+inline const mystring &ex_to_mystring(const ex &e)
+@{
+    return static_cast<const mystring &>(*e.bp);
+@}
+@end example
+
+that let's you get at the object inside an expression (after you have verified
+that the type is correct) so you can call member functions that are specific
+to the class.
+
 That's it. May the source be with you!
 
 
-@node A Comparison With Other CAS, Advantages, Symbolic functions, Top
+@node A Comparison With Other CAS, Advantages, Adding classes, Top
 @c    node-name, next, previous, up
 @chapter A Comparison With Other CAS
 @cindex advocacy
@@ -2611,6 +3689,7 @@ as you try to change the second.  Consider the simple sequence of code:
 
 @example
 #include <ginac/ginac.h>
+using namespace std;
 using namespace GiNaC;
 
 int main()
@@ -2642,6 +3721,7 @@ can be:
 
 @example
 #include <ginac/ginac.h>
+using namespace std;
 using namespace GiNaC;
 
 int main()
@@ -2935,13 +4015,13 @@ and the @samp{AM_PATH_GINAC} macro. The program used here is @file{simple.cpp}:
 
 @example
 #include <ginac/ginac.h>
-using namespace GiNaC;
 
 int main(void)
 @{
-    symbol x("x");
-    ex a = sin(x); 
-    cout << "Derivative of " << a << " is " << a.diff(x) << endl;
+    GiNaC::symbol x("x");
+    GiNaC::ex a = GiNaC::sin(x);
+    std::cout << "Derivative of " << a 
+              << " is " << a.diff(x) << std::endl;
     return 0;
 @}
 @end example
@@ -2961,9 +4041,9 @@ AC_PROG_CXX
 AC_PROG_INSTALL
 AC_LANG_CPLUSPLUS
 
-AM_PATH_GINAC(0.4.0, [
+AM_PATH_GINAC(0.7.0, [
   LIBS="$LIBS $GINACLIB_LIBS"
-  CPPFLAGS="$CFLAGS $GINACLIB_CPPFLAGS"  
+  CPPFLAGS="$CPPFLAGS $GINACLIB_CPPFLAGS"  
 ], AC_MSG_ERROR([need to have GiNaC installed]))
 
 AC_OUTPUT(Makefile)
@@ -2972,13 +4052,10 @@ AC_OUTPUT(Makefile)
 The only command in this which is not standard for automake
 is the @samp{AM_PATH_GINAC} macro.
 
-That command does the following:
-
-@display
-If a GiNaC version greater than 0.4.0 is found, adds @env{$GINACLIB_LIBS} to 
-@env{$LIBS} and @env{$GINACLIB_CPPFLAGS} to @env{$CPPFLAGS}. Otherwise, dies
-with the error message `need to have GiNaC installed'
-@end display
+That command does the following: If a GiNaC version greater or equal
+than 0.7.0 is found, then it adds @env{$GINACLIB_LIBS} to @env{$LIBS}
+and @env{$GINACLIB_CPPFLAGS} to @env{$CPPFLAGS}. Otherwise, it dies with
+the error message `need to have GiNaC installed'
 
 And the @file{Makefile.am}, which will be used to build the Makefile.
 
@@ -3059,4 +4136,3 @@ Academic Press, London
 @printindex cp
 
 @bye
-