]> www.ginac.de Git - ginac.git/blob - doc/tutorial/ginac.texi
b7ac73ad25fcf804fa632ea3eef741a7612ef067
[ginac.git] / doc / tutorial / ginac.texi
1 \input texinfo  @c -*-texinfo-*-
2 @c %**start of header
3 @setfilename ginac.info
4 @settitle GiNaC, an open framework for symbolic computation within the C++ programming language
5 @setchapternewpage off
6 @afourpaper
7 @c %**end of header
8
9 @include version.texi
10
11 @direntry
12 * ginac: (ginac).                   C++ library for symbolic computation.
13 @end direntry
14
15 @ifinfo
16 This file documents GiNaC @value{VERSION}, an open framework for symbolic
17 computation within the C++ programming language.
18
19 Copyright (C) 1999 Johannes Gutenberg University Mainz, Germany
20
21 Permission is granted to make and distribute verbatim copies of
22 this manual provided the copyright notice and this permission notice
23 are preserved on all copies.
24
25 @ignore
26 Permission is granted to process this file through TeX and print the
27 results, provided the printed document carries copying permission
28 notice identical to this one except for the removal of this paragraph
29
30 @end ignore
31 Permission is granted to copy and distribute modified versions of this
32 manual under the conditions for verbatim copying, provided that the entire
33 resulting derived work is distributed under the terms of a permission
34 notice identical to this one.
35 @end ifinfo
36
37 @titlepage
38 @title GiNaC @value{VERSION}
39 @subtitle An open framework for symbolic computation within the C++ programming language
40 @subtitle @value{UPDATED}
41 @author The GiNaC Group:
42 @author Christian Bauer, Alexander Frink, Richard B. Kreckel
43
44 @page
45 @vskip 0pt plus 1filll
46 Copyright @copyright{} 1999 Johannes Gutenberg University Mainz, Germany
47 @sp 2
48 Permission is granted to make and distribute verbatim copies of
49 this manual provided the copyright notice and this permission notice
50 are preserved on all copies.
51
52 Permission is granted to copy and distribute modified versions of this
53 manual under the conditions for verbatim copying, provided that the entire
54 resulting derived work is distributed under the terms of a permission
55 notice identical to this one.
56 @end titlepage
57
58 @page
59 @contents
60
61 @page
62
63
64 @node Top, Introduction, (dir), (dir)
65 @c    node-name, next, previous, up
66 @top GiNaC
67
68 This file documents GiNaC @value{VERSION}, an open framework for symbolic
69 computation within the C++ programming language.
70
71 @menu
72 * Introduction::                 GiNaC's purpose.
73 * A Tour of GiNaC::              A quick tour of the library.
74 * Installation::                 How to install the package.
75 * Basic Concepts::               Description of fundamental classes.
76 * Important Algorithms::         Algorithms for symbolic manipulations.
77 * Extending GiNaC::              How to extend the library.
78 * A Comparison With Other CAS::  Compares GiNaC to traditional CAS.
79 * Bibliography::
80 * Concept Index::
81 @end menu
82
83
84 @node Introduction, A Tour of GiNaC, Top, Top
85 @c    node-name, next, previous, up
86 @chapter Introduction
87
88 The motivation behind GiNaC derives from the observation that most
89 present day computer algebra systems (CAS) are linguistically and
90 semantically impoverished.  It is an attempt to overcome the current
91 situation by extending a well established and standardized computer
92 language (C++) by some fundamental symbolic capabilities, thus
93 allowing for integrated systems that embed symbolic manipulations
94 together with more established areas of computer science (like
95 computation-intense numeric applications, graphical interfaces, etc.)
96 under one roof.
97
98 This tutorial is intended for the novice user who is new to
99 GiNaC but already has some background in C++ programming.  However,
100 since a hand-made documentation like this one is difficult to keep in
101 sync with the development the actual documentation is inside the
102 sources in the form of comments.  That documentation may be parsed by
103 one of the many Javadoc-like documentation systems.  If you fail at
104 generating it you may access it from
105 @uref{http://www.ginac.de/reference/, the GiNaC home page}.
106 It is an invaluable resource not only for the advanced user who wishes
107 to extend the system (or chase bugs) but for everybody who wants to
108 comprehend the inner workings of GiNaC.  This little tutorial on the
109 other hand only covers the basic things that are unlikely to change in
110 the near future.
111
112 @section License
113 The GiNaC framework for symbolic computation within the C++ programming
114 language is Copyright @copyright{} 1999 Johannes Gutenberg University Mainz,
115 Germany.
116
117 This program is free software; you can redistribute it and/or
118 modify it under the terms of the GNU General Public License as
119 published by the Free Software Foundation; either version 2 of the
120 License, or (at your option) any later version.
121
122 This program is distributed in the hope that it will be useful, but
123 WITHOUT ANY WARRANTY; without even the implied warranty of
124 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
125 General Public License for more details.
126
127 You should have received a copy of the GNU General Public License
128 along with this program; see the file COPYING.  If not, write to the
129 Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
130 MA 02111-1307, USA.
131
132
133 @node A Tour of GiNaC, How to use it from within C++, Introduction, Top
134 @c    node-name, next, previous, up
135 @chapter A Tour of GiNaC
136
137 This quick tour of GiNaC wants to arise your interest in the
138 subsequent chapters by showing off a bit.  Please excuse us if it
139 leaves many open questions.
140
141 @menu
142 * How to use it from within C++::  Two simple examples.
143 * What it can do for you::         A Tour of GiNaC's features.
144 @end menu
145
146
147 @node How to use it from within C++, What it can do for you, A Tour of GiNaC, A Tour of GiNaC
148 @c    node-name, next, previous, up
149 @section How to use it from within C++
150
151 The GiNaC open framework for symbolic computation within the C++ programming
152 language does not try to define a language of it's own as conventional
153 CAS do.  Instead, it extends the capabilities of C++ by symbolic
154 manipulations.  Here is how to generate and print a simple (and
155 pointless) bivariate polynomial with some large coefficients:
156
157 @subheading My first GiNaC program (a bivariate polynomial)
158 @example
159 #include <ginac/ginac.h>
160 using namespace GiNaC;
161
162 int main()
163 @{
164     symbol x("x"), y("y");
165     ex poly;
166
167     for (int i=0; i<3; ++i)
168         poly += factorial(i+16)*pow(x,i)*pow(y,2-i);
169
170     cout << poly << endl;
171     return 0;
172 @}
173 @end example
174
175 Assuming the file is called @file{hello.cc}, on our system we can compile
176 and run it like this:
177
178 @example
179 $ c++ hello.cc -o hello -lcln -lginac
180 $ ./hello
181 355687428096000*x*y+20922789888000*y^2+6402373705728000*x^2
182 @end example
183
184 Next, there is a more meaningful C++ program that calls a function which
185 generates Hermite polynomials in a specified free variable.
186
187 @subheading My second GiNaC program (Hermite polynomials)
188 @example
189 #include <ginac/ginac.h>
190 using namespace GiNaC;
191
192 ex HermitePoly(symbol x, int deg)
193 @{
194     ex HKer=exp(-pow(x,2));
195     // uses the identity H_n(x) == (-1)^n exp(x^2) (d/dx)^n exp(-x^2) 
196     return normal(pow(-1,deg) * diff(HKer, x, deg) / HKer);
197 @}
198
199 int main()
200 @{
201     symbol z("z");
202
203     for (int i=0; i<6; ++i)
204         cout << "H_" << i << "(z) == " << HermitePoly(z,i) << endl;
205
206     return 0;
207 @}
208 @end example
209
210 When run, this will type out
211
212 @example
213 H_0(z) == 1
214 H_1(z) == 2*z
215 H_2(z) == 4*z^2-2
216 H_3(z) == -12*z+8*z^3
217 H_4(z) == -48*z^2+16*z^4+12
218 H_5(z) == 120*z-160*z^3+32*z^5
219 @end example
220
221 This method of generating the coefficients is of course far from
222 optimal for production purposes.
223
224 In order to show some more examples of what GiNaC can do we
225 will now use @command{ginsh}, a simple GiNaC interactive
226 shell that provides a convenient window into GiNaC's capabilities.
227
228
229 @node What it can do for you, Installation, How to use it from within C++, A Tour of GiNaC
230 @c    node-name, next, previous, up
231 @section What it can do for you
232
233 After invoking @command{ginsh} one can test and experiment with GiNaC's
234 features much like in other Computer Algebra Systems except that it does
235 not provide programming constructs like loops or conditionals.  For a
236 concise description of the @command{ginsh} syntax we refer to its
237 accompanied man page. Suffice to say that assignments and comparisons in
238 @command{ginsh} are written as they are in C, i.e. @code{=} assigns and
239 @code{==} compares.
240
241 It can manipulate arbitrary precision integers in a very fast
242 way.  Rational numbers are automatically converted to fractions of
243 coprime integers:
244
245 @example
246 > x=3^150;
247 369988485035126972924700782451696644186473100389722973815184405301748249
248 > y=3^149;
249 123329495011708990974900260817232214728824366796574324605061468433916083
250 > x/y;
251 3
252 > y/x;
253 1/3
254 @end example
255
256 All numbers occuring in GiNaC's expressions can be converted into floating
257 point numbers with the @code{evalf} method, to arbitrary accuracy:
258
259 @example
260 > evalf(1/7);
261 0.14285714285714285714
262 > Digits=150;
263 150
264 > evalf(1/7);
265 0.1428571428571428571428571428571428571428571428571428571428571428571428
266 5714285714285714285714285714285714285
267 @end example
268
269 Exact numbers other than rationals that can be manipulated in GiNaC
270 include predefined constants like Archimedes' @code{Pi}.  They can both
271 be used in symbolic manipulations (as an exact number) as well as in
272 numeric expressions (as an inexact number):
273
274 @example
275 > a=Pi^2+x;
276 x+Pi^2
277 > evalf(a);
278 x+9.869604401089358619L0
279 > x=2;
280 2
281 > evalf(a);
282 11.869604401089358619L0
283 @end example
284
285 Built-in functions evaluate immediately to exact numbers if
286 this is possible.  Conversions that can be safely performed are done
287 immediately; conversions that are not generally valid are not done:
288
289 @example
290 > cos(42*Pi);
291 1
292 > cos(acos(x));
293 x
294 > acos(cos(x));
295 acos(cos(x))
296 @end example
297
298 (Note that converting the last input to @code{x} would allow one to
299 conclude that @code{42*Pi} is equal to @code{0}.)
300
301 Linear equation systems can be solved along with basic linear
302 algebra manipulations over symbolic expressions.  In C++ GiNaC offers
303 a matrix class for this purpose but we can see what it can do using
304 @command{ginsh}'s notation of double brackets to type them in:
305
306 @example
307 > lsolve(a+x*y==z,x);
308 y^(-1)*(z-a);
309 lsolve([3*x+5*y == 7, -2*x+10*y == -5], [x, y]);
310 [x==19/8,y==-1/40]
311 > M = [[ [[1, 3]], [[-3, 2]] ]];
312 [[ [[1,3]], [[-3,2]] ]]
313 > determinant(M);
314 11
315 > charpoly(M,lambda);
316 lambda^2-3*lambda+11
317 @end example
318
319 Multivariate polynomials and rational functions may be expanded,
320 collected and normalized (i.e. converted to a ratio of two coprime 
321 polynomials):
322
323 @example
324 > a = x^4 + 2*x^2*y^2 + 4*x^3*y + 12*x*y^3 - 3*y^4;
325 -3*y^4+x^4+12*x*y^3+2*x^2*y^2+4*x^3*y
326 > b = x^2 + 4*x*y - y^2;
327 -y^2+x^2+4*x*y
328 > expand(a*b);
329 3*y^6+x^6-24*x*y^5+43*x^2*y^4+16*x^3*y^3+17*x^4*y^2+8*x^5*y
330 > collect(a*b,x);
331 3*y^6+48*x*y^4+2*x^2*y^2+x^4*(-y^2+x^2+4*x*y)+4*x^3*y*(-y^2+x^2+4*x*y)
332 > normal(a/b);
333 3*y^2+x^2
334 @end example
335
336 You can differentiate functions and expand them as Taylor or
337 Laurent series (the third argument of series is the evaluation point,
338 the fourth defines the order):
339
340 @example
341 > diff(tan(x),x);
342 tan(x)^2+1
343 > series(sin(x),x,0,4);
344 x-1/6*x^3+Order(x^4)
345 > series(1/tan(x),x,0,4);
346 x^(-1)-1/3*x+Order(x^2)
347 @end example
348
349 If you ever wanted to convert units in C or C++ and found this
350 is cumbersome, here is the solution.  Symbolic types can always be
351 used as tags for different types of objects.  Converting from wrong
352 units to the metric system is therefore easy:
353
354 @example
355 > in=.0254*m;
356 0.0254*m
357 > lb=.45359237*kg;
358 0.45359237*kg
359 > 200*lb/in^2;
360 140613.91592783185568*kg*m^(-2)
361 @end example
362
363
364 @node Installation, Prerequisites, What it can do for you, Top
365 @c    node-name, next, previous, up
366 @chapter Installation
367
368 GiNaC's installation follows the spirit of most GNU software. It is
369 easily installed on your system by three steps: configuration, build,
370 installation.
371
372 @menu
373 * Prerequisites::                Packages upon which GiNaC depends.
374 * Configuration::                How to configure GiNaC.
375 * Building GiNaC::               How to compile GiNaC.
376 * Installing GiNaC::             How to install GiNaC on your system.
377 @end menu
378
379
380 @node Prerequisites, Configuration, Installation, Installation
381 @c    node-name, next, previous, up
382 @section Prerequisites
383
384 In order to install GiNaC on your system, some prerequisites need
385 to be met.  First of all, you need to have a C++-compiler adhering to
386 the ANSI-standard @cite{ISO/IEC 14882:1998(E)}.  We used @acronym{GCC} for
387 development so if you have a different compiler you are on your own.
388 For the configuration to succeed you need a Posix compliant shell
389 installed in @file{/bin/sh}, GNU @command{bash} is fine.  Perl is needed
390 by the built process as well, since some of the source files are automatically
391 generated by Perl scripts.  Last but not least, Bruno Haible's library
392 @acronym{CLN} is extensively used and needs to be installed on your system.
393 Please get it from @uref{ftp://ftp.santafe.edu/pub/gnu/} or from
394 @uref{ftp://ftp.ilog.fr/pub/Users/haible/gnu/, Bruno Haible's FTP site}
395 (it is covered by GPL) and install it prior to trying to install GiNaC.
396 The configure script checks if it can find it and if it cannot
397 it will refuse to continue.
398
399
400 @node Configuration, Building GiNaC, Prerequisites, Installation
401 @c    node-name, next, previous, up
402 @section Configuration
403
404 To configure GiNaC means to prepare the source distribution for
405 building.  It is done via a shell script called @command{configure}
406 that is shipped with the sources. (Actually, this script is by itself
407 created with GNU Autoconf from the files @file{configure.in} and
408 @file{aclocal.m4}.)  Since a configure script generated by
409 GNU Autoconf never prompts, all customization must be done either via
410 command line parameters or environment variables.  It accepts a list
411 of parameters, the complete set of which can be listed by calling it
412 with the @option{--help} option.  The most important ones will be
413 shortly described in what follows:
414
415 @itemize @bullet
416
417 @item
418 @option{--disable-shared}: When given, this option switches off the
419 build of a shared library, i.e. a @file{.so} file.  This may be convenient
420 when developing because it considerably speeds up compilation.
421
422 @item
423 @option{--prefix=@var{PREFIX}}: The directory where the compiled library
424 and headers are installed. It defaults to @file{/usr/local} which means
425 that the library is installed in the directory @file{/usr/local/lib},
426 the header files in @file{/usr/local/include/GiNaC} and the documentation
427 (like this one) into @file{/usr/local/share/doc/GiNaC}.
428
429 @item
430 @option{--libdir=@var{LIBDIR}}: Use this option in case you want to have
431 the library installed in some other directory than
432 @file{@var{PREFIX}/lib/}.
433
434 @item
435 @option{--includedir=@var{INCLUDEDIR}}: Use this option in case you want
436 to have the header files installed in some other directory than
437 @file{@var{PREFIX}/include/ginac/}. For instance, if you specify
438 @option{--includedir=/usr/include} you will end up with the header files
439 sitting in the directory @file{/usr/include/ginac/}. Note that the
440 subdirectory @file{GiNaC} is enforced by this process in order to
441 keep the header files separated from others.  This avoids some
442 clashes and allows for an easier deinstallation of GiNaC. This ought
443 to be considered A Good Thing (tm).
444
445 @item
446 @option{--datadir=@var{DATADIR}}: This option may be given in case you
447 want to have the documentation installed in some other directory than
448 @file{@var{PREFIX}/share/doc/GiNaC/}.
449
450 @end itemize
451
452 In addition, you may specify some environment variables.
453 @env{CXX} holds the path and the name of the C++ compiler
454 in case you want to override the default in your path.  (The
455 @command{configure} script searches your path for @command{c++},
456 @command{g++}, @command{gcc}, @command{CC}, @command{cxx}
457 and @command{cc++} in that order.)  It may be very useful to
458 define some compiler flags with the @env{CXXFLAGS} environment
459 variable, like optimization, debugging information and warning
460 levels.  If omitted, it defaults to @option{-g -O2}.
461
462 The whole process is illustrated in the following two
463 examples. (Substitute @command{setenv @var{VARIABLE} @var{value}} for
464 @command{export @var{VARIABLE}=@var{value}} if the Berkeley C shell is
465 your login shell.)
466
467 @subheading Sample sessions of how to call the configure script
468
469 Simple configuration for a site-wide GiNaC library assuming everything
470 is in default paths:
471
472 @example
473 $ export CXXFLAGS="-Wall -O2"
474 $ ./configure
475 @end example
476
477 Configuration for a private static GiNaC library with several components
478 sitting in custom places (site-wide @acronym{GCC} and private @acronym{CLN}),
479 the compiler pursuaded to be picky and full assertions switched on:
480
481 @example
482 $ export CXX=/usr/local/gnu/bin/c++
483 $ export CPPFLAGS="$(CPPFLAGS) -I$(HOME)/include"
484 $ export CXXFLAGS="$(CXXFLAGS) -DDO_GINAC_ASSERT -ggdb -Wall -ansi -pedantic -O2"
485 $ export LDFLAGS="$(LDFLAGS) -L$(HOME)/lib"
486 $ ./configure --disable-shared --prefix=$(HOME)
487 @end example
488
489
490 @node Building GiNaC, Installing GiNaC, Configuration, Installation
491 @c    node-name, next, previous, up
492 @section Building GiNaC
493
494 After proper configuration you should just build the whole
495 library by typing
496
497 @example
498 $ make
499 @end example
500
501 at the command prompt and go for a cup of coffee.
502
503 Just to make sure GiNaC works properly you may run a simple test
504 suite by typing
505
506 @example
507 $ make check
508 @end example
509
510 This will compile some sample programs, run them and compare the output
511 to reference output. Each of the checks should return a message @samp{passed}
512 together with the CPU time used for that particular test.  If it does
513 not, something went wrong.  This is mostly intended to be a QA-check
514 if something was broken during the development, but not a sanity check
515 of your system.  Another intent is to allow people to fiddle around
516 with optimization.  If @acronym{CLN} was installed all right
517 this step is unlikely to return any errors.
518
519
520 @node Installing GiNaC, Basic Concepts, Building GiNaC, Installation
521 @c    node-name, next, previous, up
522 @section Installing GiNaC
523
524 To install GiNaC on your system, simply type
525
526 @example
527 $ make install
528 @end example
529
530 As described in the section about configuration
531 the files will be installed in the following directories (the
532 directories will be created if they don't already exist):
533
534 @itemize @bullet
535
536 @item
537 @file{libginac.a} will go into @file{@var{PREFIX}/lib/} (or
538 @file{@var{LIBDIR}}) which defaults to @file{/usr/local/lib/}.
539 So will @file{libginac.so} unless the configure script was
540 given the option @option{--disable-shared}.  The proper symlinks
541 will be established as well.
542
543 @item
544 All the header files will be installed into @file{@var{PREFIX}/include/ginac/}
545 (or @file{@var{INCLUDEDIR}/ginac/}, if specified).
546
547 @item
548 All documentation (HTML and Postscript) will be stuffed into
549 @file{@var{PREFIX}/share/doc/GiNaC/} (or @file{@var{DATADIR}/doc/GiNaC/}, if
550 specified).
551
552 @end itemize
553
554 Just for the record we will list some other useful make targets:
555 @command{make clean} deletes all files generated by @command{make},
556 i.e. all the object files.  In addition @command{make distclean}
557 removes all files generated by the configuration.  And finally
558 @command{make uninstall} removes the installed library and header
559 files@footnote{Uninstallation does not work after you have called
560 @command{make distclean} since the @file{Makefile} is itself generated
561 by the configuration from @file{Makefile.in} and hence deleted by
562 @command{make distclean}.  There are two obvious ways out of this
563 dilemma.  First, you can run the configuration again with the same
564 @var{PREFIX} thus creating a @file{Makefile} with a working
565 @samp{uninstall} target.  Second, you can do it by hand since you
566 now know where all the files went during installation.}.
567
568 @node Basic Concepts, Expressions, Installing GiNaC, Top
569 @c    node-name, next, previous, up
570 @chapter Basic Concepts
571
572 This chapter will describe the different fundamental objects
573 that can be handled by GiNaC.  But before doing so, it is worthwhile
574 introducing you to the more commonly used class of expressions,
575 representing a flexible meta-class for storing all mathematical
576 objects.
577
578 @menu
579 * Expressions::                  The fundamental GiNaC class.
580 * The Class Hierarchy::          Overview of GiNaC's classes.
581 * Symbols::                      Symbolic objects.
582 * Numbers::                      Numerical objects.
583 * Constants::                    Pre-defined constants.
584 * Fundamental operations::       The power, add and mul classes
585 * Built-in functions::           Mathematical functions.
586 @end menu
587
588
589 @node Expressions, The Class Hierarchy, Basic Concepts, Basic Concepts
590 @c    node-name, next, previous, up
591 @section Expressions
592
593 The most common class of objects a user deals with is the
594 expression @code{ex}, representing a mathematical object
595 like a variable, number, function, sum, product, etc...  Expressions
596 may be put together to form new expressions, passed as arguments to
597 functions, and so on.  Here is a little collection of valid
598 expressions:
599
600 @subheading Examples of expressions
601 @example
602 ex MyEx1 = 5;                       // simple number
603 ex MyEx2 = x + 2*y;                 // polynomial in x and y
604 ex MyEx3 = (x + 1)/(x - 1);         // rational expression
605 ex MyEx4 = sin(x + 2*y) + 3*z + 41; // containing a function
606 ex MyEx5 = MyEx4 + 1;               // similar to above
607 @end example
608
609 Before describing the more fundamental objects that form the building
610 blocks of expressions we'll have a quick look under the hood by
611 describing how expressions are internally managed.
612
613 @unnumberedsubsec Digression: Expressions are reference counted
614
615 An expression is extremely light-weight since internally it
616 works like a handle to the actual representation and really holds
617 nothing more than a pointer to some other object. What this means in
618 practice is that whenever you create two @code{ex} and set
619 the second equal to the first no copying process is involved. Instead,
620 the copying takes place as soon as you try to change the second.
621 Consider the simple sequence of code:
622
623 @subheading Simple copy-on-write semantics
624 @example
625 #include <ginac/ginac.h>
626 using namespace GiNaC;
627
628 int main()
629 @{
630     symbol x("x"), y("y"), z("z");
631     ex e1, e2;
632
633     e1 = sin(x + 2*y) + 3*z + 41;
634     e2 = e1;                // e2 points to same object as e1
635     cout << e2 << endl;     // prints sin(x+2*y)+3*z+41
636     e2 += 1;                // e2 is copied into a new object
637     cout << e2 << endl;     // prints sin(x+2*y)+3*z+42
638     // ...
639 @}
640 @end example
641
642 The line @code{e2 = e1;} creates a second expression pointing to the
643 object held already by @code{e1}.  The time involved for this operation
644 is therefore constant, no matter how large @code{e1} was.  Actual
645 copying, however, must take place in the line @code{e2 += 1;} because
646 @code{e1} and @code{e2} are not handles for the same object any more.
647 This concept is called @dfn{copy-on-write semantics}.  It increases
648 performance considerably whenever one object occurs multiple times and
649 represents a simple garbage collection scheme because when an @code{ex}
650 runs out of scope its destructor checks whether other expressions handle
651 the object it points to too and deletes the object from memory if that
652 turns out not to be the case.  A slightly less trivial example of
653 differentiation using the chain-rule should make clear how powerful this
654 can be.
655
656 @subheading Advanced copy-on-write semantics
657 @example
658 #include <ginac/ginac.h>
659 using namespace GiNaC;
660
661 int main()
662 @{
663     symbol x("x"), y("y");
664
665     ex e1 = x + 3*y;
666     ex e2 = pow(e1, 3);
667     ex e3 = diff(sin(e2), x);   // first derivative of sin(e2) by x
668     cout << e1 << endl          // prints x+3*y
669          << e2 << endl          // prints (x+3*y)^3
670          << e3 << endl;         // prints 3*(x+3*y)^2*cos((x+3*y)^3)
671     // ...
672 @}
673 @end example
674
675 Here, @code{e1} will actually be referenced three times while @code{e2}
676 will be referenced two times.  When the power of an expression is built,
677 that expression needs not be copied.  Likewise, since the derivative of a
678 power of an expression can be easily expressed in terms of that expression,
679 no copying of @code{e1} is involved when @code{e3} is constructed.  So,
680 when @code{e3} is constructed it will print as
681 @code{3*(x+3*y)^2*cos((x+3*y)^3)} but the argument of @code{cos()} only
682 holds a reference to @code{e2} and the factor in front is just @code{3*e1^2}.
683
684 As a user of GiNaC, you cannot see this mechanism of
685 copy-on-write semantics.  When you insert an expression into a second
686 expression, the result behaves exactly as if the contents of the first
687 expression were inserted.  But it may be useful to remember that this
688 is not what happens.  Knowing this will enable you to write much more
689 efficient code.  If you still have an uncertain feeling with
690 copy-on-write semantics, we recommend you have a look at the
691 @uref{http://www.cerfnet.com/~mpcline/c++-faq-lite/, C++-FAQ lite} by
692 Marshall Cline.  Chapter 16 covers this issue and presents an
693 implementation which is pretty close to the one in GiNaC.
694
695 So much for expressions.  But what exactly are these expressions
696 handles of?  This will be answered in the following sections.
697
698
699 @node The Class Hierarchy, Symbols, Expressions, Basic Concepts
700 @c    node-name, next, previous, up
701 @section The Class Hierarchy
702
703 GiNaC's class hierarchy consists of several classes representing
704 mathematical objects, all of which (except for @code{ex}
705 and some helpers) are internally derived from one abstract base class
706 called @code{basic}.  You do not have to deal with objects
707 of class @code{basic}, instead you'll be dealing with
708 symbols and functions of symbols.  You'll soon learn in this chapter
709 how many of the functions on symbols are really classes.  This is
710 because simple symbolic arithmetic is not supported by languages like
711 C++ so in a certain way GiNaC has to implement its own arithmetic.
712
713 To give an idea about what kinds of symbolic composits may be
714 built we have a look at the most important classes in the class
715 hierarchy.  The dashed line symbolizes a "points to" or "handles"
716 relationship while the solid lines stand for "inherits from"
717 relationships.
718
719 @subheading The GiNaC class hierarchy
720 @image{classhierarchy}
721
722 Some of the classes shown here (the ones sitting in white boxes) are
723 abstract base classes that are of no interest at all for the user.
724 They are used internally in order to avoid code duplication if
725 two or more classes derived from them share certain features.  An
726 example would be @code{expairseq}, which is a container
727 for a sequence of pairs each consisting of one expression and a number
728 (@code{numeric}).  What @emph{is} visible to the user are the derived
729 classes @code{add} and @code{mul}, representing sums of terms and products,
730 respectively.  We'll come back later to some more details about these
731 two classes and motivate the use of pairs in sums and products here.
732
733 @unnumberedsubsec Digression: Internal representation of products and sums
734
735 Although it should be completely transparent for the user of
736 GiNaC a short discussion of this topic helps to understand the sources
737 and also explain performance to a large degree.  Consider the symbolic
738 expression @math{2*d^3*(4*a+5*b-3)} which could naively be represented
739 by a tree of linear containers for addition and multiplication, one
740 container for exponentiation with base and exponent and some atomic
741 leaves of symbols and numbers in this fashion:
742
743 @subheading Naive internal representation-tree for @math{2*d^3*(4*a+5*b-3)}
744 @image{repnaive}
745
746 However, doing so results in a rather deeply nested tree which will
747 quickly become inefficient to manipulate.  If we represent the sum
748 instead as a sequence of terms, each having a purely numeric
749 multiplicative coefficient and the multiplication as a sequence of
750 terms, each having a numeric exponent, the tree becomes much more
751 flat.
752
753 @subheading Pair-wise internal representation-tree for @math{2*d^3*(4*a+5*b-3)}
754 @image{reppair}
755
756 The number @code{3} above the symbol @code{d} shows that @code{mul}
757 objects are treated similarly where the coefficients are interpreted as
758 @emph{exponents} now.  Addition of sums of terms or multiplication of
759 products with numerical exponents can be coded to be very efficient with
760 such a pair-representation.  Internally, this handling is done by many CAS in
761 this way.  It typically speeds up manipulations by an order of
762 magnitude.  The overall multiplicative factor @code{2} and
763 the additive term @code{-3} look somewhat cumbersome in
764 this representation, however, since they are still carrying a trivial
765 exponent and multiplicative factor @code{1} respectively.
766 Within GiNaC, this is avoided by adding a field that carries overall
767 numeric coefficient.
768
769 @subheading Realistic picture of GiNaC's representation-tree for @math{2*d^3*(4*a+5*b-3)}
770 @image{repreal}
771
772 This also allows for a better handling of numeric radicals, since
773 @code{sqrt(2)} can now be carried along calculations.  Now it should be
774 clear, why both classes @code{add} and @code{mul} are derived from the
775 same abstract class: the data representation is the same, only the
776 semantics differs.  In the class hierarchy, methods for polynomial
777 expansion and the like are reimplemented for @code{add} and @code{mul},
778 but the data structure is inherited from @code{expairseq}.
779
780
781 @node Symbols, Numbers, The Class Hierarchy, Basic Concepts
782 @c    node-name, next, previous, up
783 @section Symbols
784
785 Symbols are for symbolic manipulation what atoms are for
786 chemistry.  You can declare objects of class @code{symbol}
787 as any other object simply by saying @code{symbol x,y;}.
788 There is, however, a catch in here having to do with the fact that C++
789 is a compiled language.  The information about the symbol's name is
790 thrown away by the compiler but at a later stage you may want to print
791 expressions holding your symbols.  In order to avoid confusion GiNaC's
792 symbols are able to know their own name.  This is accomplished by
793 declaring its name for output at construction time in the fashion
794 @code{symbol x("x");}.  If you declare a symbol using the
795 default constructor (i.e. without string argument) the system will
796 deal out a unique name.  That name may not be suitable for printing
797 but for internal routines when no output is desired it is often
798 enough.  We'll come across examples of such symbols later in this
799 tutorial.
800
801 This implies that the strings passed to symbols at construction
802 time may not be used for comparing two of them.  It is perfectly
803 legitimate to write @code{symbol x("x"),y("x");} but it is
804 likely to lead into trouble.  Here, @code{x} and
805 @code{y} are different symbols and statements like
806 @code{x-y} will not be simplified to zero although the
807 output @code{x-x} looks funny.  Such output may also occur
808 when there are two different symbols in two scopes, for instance when
809 you call a function that declares a symbol with a name already
810 existent in a symbol in the calling function.  Again, comparing them
811 (using @code{operator==} for instance) will always reveal
812 their difference.  Watch out, please.
813
814 Although symbols can be assigned expressions for internal
815 reasons, you should not do it (and we are not going to tell you how it
816 is done).  If you want to replace a symbol with something else in an
817 expression, you can use the expression's @code{.subs()}
818 method.
819
820
821 @node Numbers, Constants, Symbols, Basic Concepts
822 @c    node-name, next, previous, up
823 @section Numbers
824
825 For storing numerical things, GiNaC uses Bruno Haible's library
826 @acronym{CLN}.  The classes therein serve as foundation
827 classes for GiNaC.  @acronym{CLN} stands for Class Library
828 for Numbers or alternatively for Common Lisp Numbers.  In order to
829 find out more about @acronym{CLN}'s internals the reader is
830 refered to the documentation of that library.  @inforef{Introduction, , cln},
831 for more information. Suffice to say that it
832 is by itself build on top of another library, the GNU Multiple
833 Precision library @acronym{GMP}, which is an extremely fast
834 library for arbitrary long integers and rationals as well as arbitrary
835 precision floating point numbers.  It is very commonly used by several
836 popular cryptographic applications.  @acronym{CLN} extends
837 @acronym{GMP} by several useful things: First, it introduces
838 the complex number field over either reals (i.e. floating point
839 numbers with arbitrary precision) or rationals.  Second, it
840 automatically converts rationals to integers if the denominator is
841 unity and complex numbers to real numbers if the imaginary part
842 vanishes and also correctly treats algebraic functions.  Third it
843 provides good implementations of state-of-the-art algorithms for all
844 trigonometric and hyperbolic functions as well as for calculation of
845 some useful constants.
846
847 The user can construct an object of class @code{numeric} in several ways.
848 The following example shows the four most important constructors: construction
849 from C-integer, construction of fractions from two integers, construction
850 from C-float and construction from a string.
851
852 @subheading Construction of numbers
853 @example
854 #include <ginac/ginac.h>
855 using namespace GiNaC;
856
857 int main()
858 @{
859     numeric two(2);                     // exact integer 2
860     numeric r(2,3);                     // exact fraction 2/3
861     numeric e(2.71828);                 // floating point number
862     numeric p("3.1415926535897932385"); // floating point number
863
864     cout << two*p << endl;  // floating point 6.283...
865     // ...
866 @}
867 @end example
868
869 Note that all those constructors are @emph{explicit} which means you are
870 not allowed to write @code{numeric two=2;}.  This is because the basic
871 objects to be handled by GiNaC are the expressions @code{ex} and we want
872 to keep things simple and wish objects like @code{pow(x,2)} to be
873 handled the same way as @code{pow(x,a)}, which means that we need to allow
874 a general @code{ex} as base and exponent.  Therefore there is an implicit
875 constructor from C-integers directly to expressions handling numerics at
876 work in most of our examples.  This design really becomes convenient when
877 one declares own functions having more than one parameter but it forbids
878 using implicit constructors because that would lead to ambiguities. 
879
880 It may be tempting to construct numbers writing @code{numeric r(3/2)}.
881 This would, however, call C's built-in operator @code{/} for integers
882 first and result in a numeric holding a plain integer 1.  @strong{Never use
883 @code{/} on integers!} Use the constructor from two integers instead, as
884 shown in the example above.  Writing @code{numeric(1)/2} may look funny but
885 works also. 
886
887 We have seen now the distinction between exact numbers and
888 floating point numbers.  Clearly, the user should never have to worry
889 about dynamically created exact numbers, since their "exactness"
890 always determines how they ought to be handled.  The situation is
891 different for floating point numbers.  Their accuracy is handled by
892 one @emph{global} variable, called @code{Digits}.  (For those readers
893 who know about Maple: it behaves very much like Maple's @code{Digits}).
894 All objects of class numeric that are constructed from then on will be
895 stored with a precision matching that number of decimal digits:
896
897 @subheading Controlling the precision of floating point numbers
898 @example
899 #include <ginac/ginac.h>
900 using namespace GiNaC;
901
902 void foo()
903 @{
904     numeric three(3.0), one(1.0);
905     numeric x = one/three;
906
907     cout << "in " << Digits << " digits:" << endl;
908     cout << x << endl;
909     cout << Pi.evalf() << endl;
910 @}
911
912 int main()
913 @{
914     foo();
915     Digits = 60;
916     foo();
917     return 0;
918 @}
919 @end example
920
921 The above example prints the following output to screen:
922
923 @example
924 in 17 digits:
925 0.333333333333333333
926 3.14159265358979324
927 in 60 digits:
928 0.333333333333333333333333333333333333333333333333333333333333333333
929 3.14159265358979323846264338327950288419716939937510582097494459231
930 @end example
931
932 It should be clear that objects of class @code{numeric} should be used
933 for constructing numbers or for doing arithmetic with them.  The objects
934 one deals with most of the time are the polymorphic expressions @code{ex}.
935
936 @subsection Tests on numbers
937
938 Once you have declared some numbers, assigned them to
939 expressions and done some arithmetic with them it is frequently
940 desired to retrieve some kind of information from them like asking
941 whether that number is integer, rational, real or complex.  For those
942 cases GiNaC provides several useful methods.  (Internally, they fall
943 back to invocations of certain CLN functions.)
944
945 As an example, let's construct some rational number, multiply it
946 with some multiple of its denominator and check what comes out:
947
948 @subheading Sample test on objects of type numeric
949 @example
950 #include <ginac/ginac.h>
951 using namespace GiNaC;
952
953 // some very important constants:
954 const numeric twentyone(21);
955 const numeric ten(10);
956 const numeric five(5);
957
958 int main()
959 @{
960     numeric answer = twentyone;
961
962     answer /= five;
963     cout << answer.is_integer() << endl;  // false, it's 21/5
964     answer *= ten;
965     cout << answer.is_integer() << endl;  // true, it's 42 now!
966     // ...
967 @}
968 @end example
969
970 Note that the variable @code{answer} is constructed here
971 as an integer by @code{numeric}'s copy constructor but in
972 an intermediate step it holds a rational number represented as integer
973 numerator and integer denominator.  When multiplied by 10, the
974 denominator becomes unity and the result is automatically converted to
975 a pure integer again.  Internally, the underlying
976 @acronym{CLN} is responsible for this behaviour and we refer
977 the reader to @acronym{CLN}'s documentation.  Suffice to say
978 that the same behaviour applies to complex numbers as well as return
979 values of certain functions.  Complex numbers are automatically
980 converted to real numbers if the imaginary part becomes zero.  The
981 full set of tests that can be applied is listed in the following
982 table.
983
984 @cartouche
985 @multitable @columnfractions .33 .66
986 @item Method @tab Returns true if@dots{}
987 @item @code{.is_zero()}
988 @tab object is equal to zero
989 @item @code{.is_positive()}
990 @tab object is not complex and greater than 0
991 @item @code{.is_integer()}
992 @tab object is a (non-complex) integer
993 @item @code{.is_pos_integer()}
994 @tab object is an integer and greater than 0
995 @item @code{.is_nonneg_integer()}
996 @tab object is an integer and greater equal 0
997 @item @code{.is_even()}
998 @tab object is an even integer
999 @item @code{.is_odd()}
1000 @tab object is an odd integer
1001 @item @code{.is_prime()}
1002 @tab object is a prime integer (probabilistic primality test)
1003 @item @code{.is_rational()}
1004 @tab object is an exact rational number (integers are rational, too, as are complex extensions like @math{2/3+7/2*I})
1005 @item @code{.is_real()}
1006 @tab object is a real integer, rational or float (i.e. is not complex)
1007 @end multitable
1008 @end cartouche
1009
1010
1011 @node Constants, Fundamental operations, Numbers, Basic Concepts
1012 @c    node-name, next, previous, up
1013 @section Constants
1014
1015 Constants behave pretty much like symbols except that they return
1016 some specific number when the method @code{.evalf()} is called.
1017
1018 The predefined known constants are:
1019
1020 @cartouche
1021 @multitable @columnfractions .14 .29 .57
1022 @item Name @tab Common Name @tab Numerical Value (35 digits)
1023 @item @code{Pi}
1024 @tab Archimedes' constant
1025 @tab 3.14159265358979323846264338327950288
1026 @item @code{Catalan}
1027 @tab Catalan's constant
1028 @tab 0.91596559417721901505460351493238411
1029 @item @code{EulerGamma}
1030 @tab Euler's (or Euler-Mascheroni) constant
1031 @tab 0.57721566490153286060651209008240243
1032 @end multitable
1033 @end cartouche
1034
1035
1036 @node Fundamental operations, Built-in functions, Constants, Basic Concepts
1037 @c    node-name, next, previous, up
1038 @section Fundamental operations: the @code{power}, @code{add} and @code{mul} classes
1039
1040 Simple polynomial expressions are written down in GiNaC pretty
1041 much like in other CAS.  The necessary operators @code{+}, @code{-},
1042 @code{*} and @code{/} have been overloaded to achieve this goal.
1043 When you run the following program, the constructor for an object of
1044 type @code{mul} is automatically called to hold the product of @code{a}
1045 and @code{b} and then the constructor for an object of type @code{add}
1046 is called to hold the sum of that @code{mul} object and the number one:
1047
1048 @subheading Construction of @code{add} and @code{mul} objects
1049 @example
1050 #include <ginac/ginac.h>
1051 using namespace GiNaC;
1052
1053 int main()
1054 @{
1055     symbol a("a"), b("b");
1056     ex MyTerm = 1+a*b;
1057     // ...
1058 @}
1059 @end example
1060
1061 For exponentiation, you have already seen the somewhat clumsy (though C-ish)
1062 statement @code{pow(x,2);} to represent @code{x} squared.  This direct
1063 construction is necessary since we cannot safely overload the constructor
1064 @code{^} in C++ to construct a @code{power} object.  If we did, it would
1065 have several counterintuitive effects:
1066
1067 @itemize @bullet
1068 @item
1069 Due to C's operator precedence, @code{2*x^2} would be parsed as @code{(2*x)^2}.
1070 @item
1071 Due to the binding of the operator @code{^}, @code{x^a^b} would result in
1072 @code{(x^a)^b}. This would be confusing since most (though not all) other CAS
1073 interpret this as @code{x^(a^b)}.
1074 @item
1075 Also, expressions involving integer exponents are very frequently used,
1076 which makes it even more dangerous to overload @code{^} since it is then
1077 hard to distinguish between the semantics as exponentiation and the one
1078 for exclusive or.  (It would be embarassing to return @code{1} where one
1079 has requested @code{2^3}.)
1080 @end itemize
1081
1082 All effects are contrary to mathematical notation and differ from the
1083 way most other CAS handle exponentiation, therefore overloading
1084 @code{^} is ruled out for GiNaC's C++ part.  The situation
1085 is different in @command{ginsh}, there the exponentiation-@code{^}
1086 exists.  (Also note, that the other frequently used exponentiation operator
1087 @code{**} does not exist at all in C++).
1088
1089 To be somewhat more precise, objects of the three classes
1090 described here, are all containers for other expressions.  An object
1091 of class @code{power} is best viewed as a container with
1092 two slots, one for the basis, one for the exponent.  All valid GiNaC
1093 expressions can be inserted.  However, basic transformations like
1094 simplifying @code{pow(pow(x,2),3)} to @code{x^6} automatically are only
1095 performed when this is mathematically possible.  If we replace the
1096 outer exponent three in the example by some symbols @code{a}, the
1097 simplification is not safe and will not be performed, since @code{a}
1098 might be @code{1/2} and @code{x} negative.
1099
1100 Objects of type @code{add} and @code{mul} are containers with an arbitrary
1101 number of slots for expressions to be inserted.  Again, simple and safe
1102 simplifications are carried out like transforming @code{3*x+4-x} to
1103 @code{2*x+4}.
1104
1105 The general rule is that when you construct such objects, GiNaC
1106 automatically creates them in canonical form, which might differ from
1107 the form you typed in your program.  This allows for rapid comparison
1108 of expressions, since after all @code{a-a} is simply zero.
1109 Note, that the canonical form is not necessarily lexicographical
1110 ordering or in any way easily guessable.  It is only guaranteed that
1111 constructing the same expression twice, either implicitly or
1112 explicitly, results in the same canonical form.
1113
1114
1115 @node Built-in functions, Important Algorithms, Fundamental operations, Basic Concepts
1116 @c    node-name, next, previous, up
1117 @section Built-in functions
1118
1119 There are quite a number of useful functions built into GiNaC.
1120 They are all objects of class @code{function}.  They
1121 accept one or more expressions as arguments and return one expression.
1122 If the arguments are not numerical, the evaluation of the functions
1123 may be halted, as it does in the next example:
1124
1125 @subheading Evaluation of built-in functions
1126 @example
1127 #include <ginac/ginac.h>
1128 using namespace GiNaC;
1129
1130 int main()
1131 @{
1132     symbol x("x"), y("y");
1133     
1134     ex foo = x+y/2;
1135     cout << "gamma(" << foo << ") -> " << gamma(foo) << endl;
1136     ex bar = foo.subs(y==1);
1137     cout << "gamma(" << bar << ") -> " << gamma(bar) << endl;
1138     ex foobar = bar.subs(x==7);
1139     cout << "gamma(" << foobar << ") -> " << gamma(foobar) << endl;
1140     // ...
1141 @}
1142 @end example
1143
1144 This program will type out two times a function and then an
1145 expression that may be really useful:
1146
1147 @example
1148 gamma(x+(1/2)*y) -> gamma(x+(1/2)*y)
1149 gamma(x+1/2) -> gamma(x+1/2)
1150 gamma(15/2) -> (135135/128)*Pi^(1/2)
1151 @end example
1152
1153 Most of these functions can be differentiated, series expanded and so
1154 on.  Read the next chapter in order to learn more about this..
1155
1156
1157 @node Important Algorithms, Polynomial Expansion, Built-in functions, Top
1158 @c    node-name, next, previous, up
1159 @chapter Important Algorithms
1160
1161 In this chapter the most important algorithms provided by GiNaC
1162 will be described.  Some of them are implemented as functions on
1163 expressions, others are implemented as methods provided by expression
1164 objects.  If they are methods, there exists a wrapper function around
1165 it, so you can alternatively call it in a functional way as shown in
1166 the simple example:
1167
1168 @subheading Methods vs. wrapper functions
1169 @example
1170 #include <ginac/ginac.h>
1171 using namespace GiNaC;
1172
1173 int main()
1174 @{
1175     ex x = numeric(1.0);
1176     
1177     cout << "As method:   " << sin(x).evalf() << endl;
1178     cout << "As function: " << evalf(sin(x)) << endl;
1179     // ...
1180 @}
1181 @end example
1182
1183 The general rule is that wherever methods accept one or more
1184 parameters (@var{arg1}, @var{arg2}, @dots{}) the order of arguments
1185 the function wrapper accepts is the same but preceded by the object
1186 to act on (@var{object}, @var{arg1}, @var{arg2}, @dots{}).  This
1187 approach is the most natural one in an OO model but it may lead to
1188 confusion for MapleV users because where they would type
1189 @code{A:=x+1; subs(x=2,A);} GiNaC would require
1190 @code{A=x+1; subs(A,x==2);} (after proper declaration of @code{A}
1191 and @code{x}).  On the other hand, since MapleV returns 3 on
1192 @code{A:=x^2+3; coeff(A,x,0);} (GiNaC:
1193 @code{A=pow(x,2)+3; coeff(A,x,0);}) it is clear that
1194 MapleV is not trying to be consistent here.  Also, users of MuPAD will
1195 in most cases feel more comfortable with GiNaC's convention.  All
1196 function wrappers are always implemented as simple inline functions
1197 which just call the corresponding method and are only provided for
1198 users uncomfortable with OO who are dead set to avoid method
1199 invocations.  Generally, a chain of function wrappers is much harder
1200 to read than a chain of methods and should therefore be avoided if
1201 possible.  On the other hand, not everything in GiNaC is a method on
1202 class @code{ex} and sometimes calling a function cannot be
1203 avoided.
1204
1205 @menu
1206 * Polynomial Expansion::
1207 * Collecting expressions::
1208 * Polynomial Arithmetic::
1209 * Symbolic Differentiation::
1210 * Series Expansion::
1211 @end menu
1212
1213
1214 @node Polynomial Expansion, Collecting expressions, Important Algorithms, Important Algorithms
1215 @c    node-name, next, previous, up
1216 @section Polynomial Expansion
1217
1218 A polynomial in one or more variables has many equivalent
1219 representations.  Some useful ones serve a specific purpose.  Consider
1220 for example the trivariate polynomial @math{4*x*y + x*z + 20*y^2 + 21*y*z + 4*z^2}.
1221 It is equivalent to the factorized polynomial @math{(x + 5*y + 4*z)*(4*y + z)}.
1222 Other representations are the recursive ones where one collects for
1223 exponents in one of the three variable.  Since the factors are
1224 themselves polynomials in the remaining two variables the procedure
1225 can be repeated.  In our expample, two possibilies would be
1226 @math{(4*y + z)*x + 20*y^2 + 21*y*z + 4*z^2} and
1227 @math{20*y^2 + (21*z + 4*x)*y + 4*z^2 + x*z}.
1228
1229 To bring an expression into expanded form, its method
1230 @code{.expand()} may be called.  In our example above,
1231 this corresponds to @math{4*x*y + x*z + 20*y^2 + 21*y*z + 4*z^2}.
1232 Again, since the canonical form in GiNaC is not easily guessable you
1233 should be prepared to see different orderings of terms in such sums!
1234
1235
1236 @node Collecting expressions, Polynomial Arithmetic, Polynomial Expansion, Important Algorithms
1237 @c    node-name, next, previous, up
1238 @section Collecting expressions
1239
1240 Another useful representation of multivariate polynomials is as
1241 a univariate polynomial in one of the variables with the coefficients
1242 being polynomials in the remaining variables.  The method
1243 @code{collect()} accomplishes this task:
1244
1245 @example
1246 #include <ginac/ginac.h>
1247 ex ex::collect(symbol const & s);
1248 @end example
1249
1250 Note that the original polynomial needs to be in expanded form in
1251 order to be able to find the coefficients properly.  The range of
1252 occuring coefficients can be checked using the two methods
1253
1254 @example
1255 #include <ginac/ginac.h>
1256 int ex::degree(symbol const & s);
1257 int ex::ldegree(symbol const & s);
1258 @end example
1259
1260 where @code{degree()} returns the highest coefficient and
1261 @code{ldegree()} the lowest one.  These two methods work
1262 also reliably on non-expanded input polynomials.  This is illustrated
1263 in the following example: 
1264
1265 @subheading Collecting expressions in multivariate polynomials
1266 @example
1267 #include <ginac/ginac.h>
1268 using namespace GiNaC;
1269
1270 int main()
1271 @{
1272     symbol x("x"), y("y");
1273     ex PolyInp = 4*pow(x,3)*y + 5*x*pow(y,2) + 3*y
1274                  - pow(x+y,2) + 2*pow(y+2,2) - 8;
1275     ex Poly = PolyInp.expand();
1276     
1277     for (int i=Poly.ldegree(x); i<=Poly.degree(x); ++i) @{
1278         cout << "The x^" << i << "-coefficient is "
1279              << Poly.coeff(x,i) << endl;
1280     @}
1281     cout << "As polynomial in y: " 
1282          << Poly.collect(y) << endl;
1283     // ...
1284 @}
1285 @end example
1286
1287 When run, it returns an output in the following fashion:
1288
1289 @example
1290 The x^0-coefficient is y^2+11*y
1291 The x^1-coefficient is 5*y^2-2*y
1292 The x^2-coefficient is -1
1293 The x^3-coefficient is 4*y
1294 As polynomial in y: -x^2+(5*x+1)*y^2+(-2*x+4*x^3+11)*y
1295 @end example
1296
1297 As always, the exact output may vary between different versions of
1298 GiNaC or even from run to run since the internal canonical ordering is
1299 not within the user's sphere of influence.
1300
1301
1302 @node Polynomial Arithmetic, Symbolic Differentiation, Collecting expressions, Important Algorithms
1303 @c    node-name, next, previous, up
1304 @section Polynomial Arithmetic
1305
1306 @subsection GCD and LCM
1307
1308 The functions for polynomial greatest common divisor and least common
1309 multiple have the synopsis:
1310
1311 @example
1312 #include <ginac/normal.h>
1313 ex gcd(const ex & a, const ex & b);
1314 ex lcm(const ex & a, const ex & b);
1315 @end example
1316
1317 The functions @code{gcd()} and @code{lcm()} accept two expressions
1318 @code{a} and @code{b} as arguments and return
1319 a new expression, their greatest common divisor or least common
1320 multiple, respectively.  If the polynomials @code{a} and
1321 @code{b} are coprime @code{gcd(a,b)} returns 1 and @code{lcm(a,b)}
1322 returns the product of @code{a} and @code{b}.
1323
1324 @subheading Polynomal GCD/LCM
1325 @example
1326 #include <ginac/ginac.h>
1327 using namespace GiNaC;
1328
1329 int main()
1330 @{
1331     symbol x("x"), y("y"), z("z");
1332     ex P_a = 4*x*y + x*z + 20*pow(y, 2) + 21*y*z + 4*pow(z, 2);
1333     ex P_b = x*y + 3*x*z + 5*pow(y, 2) + 19*y*z + 12*pow(z, 2);
1334
1335     ex P_gcd = gcd(P_a, P_b);
1336     // x + 5*y + 4*z
1337     ex P_lcm = lcm(P_a, P_b);
1338     // 4*x*y^2 + 13*y*x*z + 20*y^3 + 81*y^2*z + 67*y*z^2 + 3*x*z^2 + 12*z^3
1339     // ...
1340 @}
1341 @end example
1342
1343 @subsection The @code{normal} method
1344
1345 While in common symbolic code @code{gcd()} and @code{lcm()} are not too
1346 heavily used, simplification is called for frequently.  Therefore
1347 @code{.normal()}, which provides some basic form of simplification, has
1348 become a method of class @code{ex}, just like @code{.expand()}.  It
1349 converts a rational function into an equivalent rational function where
1350 numerator and denominator are coprime.  This means, it finds the GCD of
1351 numerator and denominator and cancels it.  If it encounters some object
1352 which does not belong to the domain of rationals (a function for
1353 instance), that object is replaced by a temporary symbol.  This means
1354 that both expressions @code{t1} and @code{t2} are indeed simplified in
1355 this little program:
1356
1357 @subheading Cancellation of polynomial GCD (with obstacles)
1358 @example
1359 #include <ginac/ginac.h>
1360 using namespace GiNaC;
1361
1362 int main()
1363 @{
1364     symbol x("x");
1365     ex t1 = (pow(x,2) + 2*x + 1)/(x + 1);
1366     ex t2 = (pow(sin(x),2) + 2*sin(x) + 1)/(sin(x) + 1);
1367     cout << "t1 is " << t1.normal() << endl;
1368     cout << "t2 is " << t2.normal() << endl;
1369     // ...
1370 @}
1371 @end example
1372
1373 Of course this works for multivariate polynomials too, so the ratio of
1374 the sample-polynomials from the section about GCD and LCM above would
1375 be normalized to @code{P_a/P_b} = @code{(4*y+z)/(y+3*z)}.
1376
1377
1378 @node Symbolic Differentiation, Series Expansion, Polynomial Arithmetic, Important Algorithms
1379 @c    node-name, next, previous, up
1380 @section Symbolic Differentiation
1381
1382 GiNaC's objects know how to differentiate themselves.  Thus, a
1383 polynomial (class @code{add}) knows that its derivative is
1384 the sum of the derivatives of all the monomials:
1385
1386 @subheading Simple polynomial differentiation
1387 @example
1388 #include <ginac/ginac.h>
1389 using namespace GiNaC;
1390
1391 int main()
1392 @{
1393     symbol x("x"), y("y"), z("z");
1394     ex P = pow(x, 5) + pow(x, 2) + y;
1395
1396     cout << P.diff(x,2) << endl;  // 20*x^3 + 2
1397     cout << P.diff(y) << endl;    // 1
1398     cout << P.diff(z) << endl;    // 0
1399     // ...
1400 @}
1401 @end example
1402
1403 If a second integer parameter @var{n} is given, the @code{diff} method
1404 returns the @var{n}th derivative.
1405
1406 If @emph{every} object and every function is told
1407 what its derivative is, all derivatives of composed objects can be
1408 calculated using the chain rule and the product rule.  Consider, for
1409 instance the expression @code{1/cosh(x)}.  Since the
1410 derivative of @code{cosh(x)} is @code{sinh(x)}
1411 and the derivative of @code{pow(x,-1)} is
1412 @code{-pow(x,-2)}, GiNaC can readily compute the
1413 composition.  It turns out that the composition is the generating
1414 function for Euler Numbers, i.e. the so called
1415 @var{n}th Euler number is the coefficient of @code{x^n/n!} in
1416 the expansion of @code{1/cosh(x)}.  We may use this identity to code a
1417 function that generates Euler numbers in just three lines:
1418
1419 @subheading Differentiation with nontrivial functions: Euler numbers
1420 @example
1421 #include <ginac/ginac.h>
1422 using namespace GiNaC;
1423
1424 ex EulerNumber(unsigned n)
1425 @{
1426     symbol x;
1427     ex generator = pow(cosh(x),-1);
1428     return generator.diff(x,n).subs(x==0);
1429 @}
1430
1431 int main()
1432 @{
1433     for (unsigned i=0; i<11; i+=2)
1434         cout << EulerNumber(i) << endl;
1435     return 0;
1436 @}
1437 @end example
1438
1439 When you run it, it produces the sequence @code{1}, @code{-1}, @code{5},
1440 @code{-61}, @code{1385}, @code{-50521}.  We increment the loop variable
1441 @code{i} by two since all odd Euler numbers vanish anyways.
1442
1443
1444 @node Series Expansion, Extending GiNaC, Symbolic Differentiation, Important Algorithms
1445 @c    node-name, next, previous, up
1446 @section Series Expansion
1447
1448 Expressions know how to expand themselves as a Taylor series or (more
1449 generally) a Laurent series.  Similar to most conventional Computer
1450 Algebra Systems, no distinction is made between those two.  There is a
1451 class of its own for storing such series as well as a class for storing
1452 the order of the series.  A sample program could read:
1453
1454 @subheading Series expansion
1455 @example
1456 #include <ginac/ginac.h>
1457 using namespace GiNaC;
1458
1459 int main()
1460 @{
1461     symbol x("x");
1462     numeric point(0);
1463     ex MyExpr1 = sin(x);
1464     ex MyExpr2 = 1/(x - pow(x, 2) - pow(x, 3));
1465     ex MyTailor, MySeries;
1466     
1467     MyTailor = MyExpr1.series(x, point, 5);
1468     cout << MyExpr1 << " == " << MyTailor
1469          << " for small " << x << endl;
1470     MySeries = MyExpr2.series(x, point, 7);
1471     cout << MyExpr2 << " == " << MySeries
1472          << " for small " << x << endl;
1473     // ...
1474 @}
1475 @end example
1476
1477 As an instructive application, let us calculate the numerical
1478 value of Archimedes' constant (for which there already exists the
1479 exbuilt-in constant @code{Pi}) using M@'echain's
1480 mysterious formula @code{Pi==16*atan(1/5)-4*atan(1/239)}.
1481 We may expand the arcus tangent around @code{0} and insert
1482 the fractions @code{1/5} and @code{1/239}.
1483 But, as we have seen, a series in GiNaC carries an order term with it.
1484 The function @code{series_to_poly()} may be used to strip 
1485 this off:
1486
1487 @subheading Series expansion using M@'echain's formula for @code{Pi}
1488 @example
1489 #include <ginac/ginac.h>
1490 using namespace GiNaC;
1491
1492 ex mechain_pi(int degr)
1493 @{
1494     symbol x;
1495     ex pi_expansion = series_to_poly(atan(x).series(x,0,degr));
1496     ex pi_approx = 16*pi_expansion.subs(x==numeric(1,5))
1497                    -4*pi_expansion.subs(x==numeric(1,239));
1498     return pi_approx;
1499 @}
1500
1501 int main()
1502 @{
1503     ex pi_frac;
1504     for (int i=2; i<12; i+=2) @{
1505         pi_frac = mechain_pi(i);
1506         cout << i << ":\t" << pi_frac << endl
1507              << "\t" << pi_frac.evalf() << endl;
1508     @}
1509     return 0;
1510 @}
1511 @end example
1512
1513 When you run this program, it will type out:
1514
1515 @example
1516 2:      3804/1195
1517         3.1832635983263598326
1518 4:      5359397032/1706489875
1519         3.1405970293260603143
1520 6:      38279241713339684/12184551018734375
1521         3.141621029325034425
1522 8:      76528487109180192540976/24359780855939418203125
1523         3.141591772182177295
1524 10:     327853873402258685803048818236/104359128170408663038552734375
1525         3.1415926824043995174
1526 @end example
1527
1528
1529 @node Extending GiNaC, What does not belong into GiNaC, Series Expansion, Top
1530 @c    node-name, next, previous, up
1531 @chapter Extending GiNaC
1532
1533 By reading so far you should have gotten a fairly good
1534 understanding of GiNaC's design-patterns.  From here on you should
1535 start reading the sources.  All we can do now is issue some
1536 recommendations how to tackle GiNaC's many loose ends in order to
1537 fulfill everybody's dreams.
1538
1539 @menu
1540 * What does not belong into GiNaC::  What to avoid.
1541 * Symbolic functions::               Implementing symbolic functions.
1542 @end menu
1543
1544
1545 @node What does not belong into GiNaC, Symbolic functions, Extending GiNaC, Extending GiNaC
1546 @c    node-name, next, previous, up
1547 @section What doesn't belong into GiNaC
1548
1549 First of all, GiNaC's name must be read literally.  It is
1550 designed to be a library for use within C++.  The tiny
1551 @command{ginsh} accompanying GiNaC makes this even more
1552 clear: it doesn't even attempt to provide a language.  There are no
1553 loops or conditional expressions in @command{ginsh}, it is
1554 merely a window into the library for the programmer to test stuff (or
1555 to show off).  Still, the design of a complete CAS with a language of
1556 its own, graphical capabilites and all this on top of GiNaC is
1557 possible and is without doubt a nice project for the future.
1558
1559 There are many built-in functions in GiNaC that do not know how
1560 to evaluate themselves numerically to a precision declared at runtime
1561 (using @code{Digits}).  Some may be evaluated at certain
1562 points, but not generally.  This ought to be fixed.  However, doing
1563 numerical computations with GiNaC's quite abstract classes is doomed
1564 to be inefficient.  For this purpose, the underlying bignum-package
1565 @acronym{CLN} is much better suited.
1566
1567
1568 @node Symbolic functions, A Comparison With Other CAS, What does not belong into GiNaC, Extending GiNaC
1569 @c    node-name, next, previous, up
1570 @section Symbolic functions
1571
1572 The easiest and most instructive way to start with is probably
1573 to implement your own function.  Objects of class
1574 @code{function} are inserted into the system via a kind of
1575 "registry".  They get a serial number that is used internally to
1576 identify them but you usually need not worry about this.  What you
1577 have to care for are functions that are called when the user invokes
1578 certain methods.  These are usual C++-functions accepting a number of
1579 @code{ex} as arguments and returning one
1580 @code{ex}.  As an example, if we have a look at a
1581 simplified implementation of the cosine trigonometric function, we
1582 first need a function that is called when one wishes to
1583 @code{eval} it.  It could look something like this:
1584
1585 @example
1586 static ex cos_eval_method(ex const & x)
1587 @{
1588     // if x%2*Pi return 1
1589     // if x%Pi return -1
1590     // if x%Pi/2 return 0
1591     // care for other cases...
1592     return cos(x).hold();
1593 @}
1594 @end example
1595
1596 The last line returns @code{cos(x)} if we don't know what
1597 else to do and stops a potential recursive evaluation by saying
1598 @code{.hold()}.  We should also implement a method for
1599 numerical evaluation and since we are lazy we sweep the problem under
1600 the rug by calling someone else's function that does so, in this case
1601 the one in class @code{numeric}:
1602
1603 @example
1604 static ex cos_evalf_method(ex const & x)
1605 @{
1606     return sin(ex_to_numeric(x));
1607 @}
1608 @end example
1609
1610 Differentiation will surely turn up and so we need to tell
1611 @code{sin} how to differentiate itself:
1612
1613 @example
1614 static ex cos_diff_method(ex const & x, unsigned diff_param)
1615 @{
1616     return cos(x);
1617 @}
1618 @end example
1619
1620 The second parameter is obligatory but uninteresting at this point.
1621 It is used for correct handling of the product rule only.  For Taylor
1622 expansion, it is enough to know how to differentiate.  But if the
1623 function you want to implement does have a pole somewhere in the
1624 complex plane, you need to write another method for Laurent expansion
1625 around that point.
1626
1627 Now that everything has been written for @code{cos},
1628 we need to tell the system about it.  This is done by a macro and we
1629 are not going to descibe how it expands, please consult your
1630 preprocessor if you are curious:
1631
1632 @example
1633 REGISTER_FUNCTION(cos, cos_eval_method, cos_evalf_method, cos_diff, NULL);
1634 @end example
1635
1636 The first argument is the function's name, the second, third and
1637 fourth bind the corresponding methods to this objects and the fifth is
1638 a slot for inserting a method for series expansion.  Also, the new
1639 function needs to be declared somewhere.  This may also be done by a
1640 convenient preprocessor macro:
1641
1642 @example
1643 DECLARE_FUNCTION_1P(cos)
1644 @end example
1645
1646 The suffix @code{_1P} stands for @emph{one parameter}.  Of course, this
1647 implementation of @code{cos} is very incomplete and lacks several safety
1648 mechanisms.  Please, have a look at the real implementation in GiNaC.
1649 (By the way: in case you are worrying about all the macros above we
1650 can assure you that functions are GiNaC's most macro-intense classes.
1651 We have done our best to avoid them where we can.)
1652
1653 That's it. May the source be with you!
1654
1655
1656 @node A Comparison With Other CAS, Bibliography, Symbolic functions, Top
1657 @c    node-name, next, previous, up
1658 @chapter A Comparison With Other CAS
1659
1660 This chapter will give you some information on how GiNaC
1661 compares to other, traditional Computer Algebra Systems, like
1662 @emph{Maple}, @emph{Mathematica} or @emph{Reduce}, where it has
1663 advantages and disadvantages over these systems.
1664
1665
1666 @heading Advantages
1667
1668 GiNaC has several advantages over traditional Computer
1669 Algebra Systems, like 
1670
1671 @itemize @bullet
1672
1673 @item
1674 familiar language: all common CAS implement their own
1675 proprietary grammar which you have to learn first (and maybe learn
1676 again when your vendor chooses to "enhance" it).  With GiNaC you
1677 can write your program in common C++, which is standardized.
1678
1679 @item
1680 structured data types: you can build up structured data
1681 types using @code{struct}s or @code{class}es
1682 together with STL features instead of using unnamed lists of lists
1683 of lists.
1684
1685 @item
1686 strongly typed: in CAS, you usually have only one kind of
1687 variables which can hold contents of an arbitrary type.  This
1688 4GL like feature is nice for novice programmers, but dangerous.
1689     
1690 @item
1691 development tools: powerful development tools exist for
1692 C++, like fancy editors (e.g. with automatic
1693 indentation and syntax highlighting), debuggers, visualization
1694 tools, documentation tools...
1695
1696 @item
1697 modularization: C++ programs can 
1698 easily be split into modules by separating interface and
1699 implementation.
1700
1701 @item
1702 price: GiNaC is distributed under the GNU Public License
1703 which means that it is free and available with source code.  And
1704 there are excellent C++-compilers for free, too.
1705     
1706 @item
1707 extendable: you can add your own classes to GiNaC, thus
1708 extending it on a very low level.  Compare this to a traditional
1709 CAS that you can usually only extend on a high level by writing in
1710 the language defined by the parser.  In particular, it turns out
1711 to be almost impossible to fix bugs in a traditional system.
1712
1713 @item
1714 seemless integration: it is somewhere between difficult
1715 and impossible to call CAS functions from within a program 
1716 written in C++ or any other programming 
1717 language and vice versa.  With GiNaC, your symbolic routines
1718 are part of your program.  You can easily call third party
1719 libraries, e.g. for numerical evaluation or graphical 
1720 interaction.  All other approaches are much more cumbersome: they
1721 range from simply ignoring the problem (i.e. @emph{Maple}) to providing a
1722 method for "embedding" the system (i.e. @emph{Yacas}).
1723
1724 @item
1725 efficiency: often large parts of a program do not need
1726 symbolic calculations at all.  Why use large integers for loop
1727 variables or arbitrary precision arithmetics where double
1728 accuracy is sufficient?  For pure symbolic applications,
1729 GiNaC is comparable in speed with other CAS.
1730
1731 @end itemize
1732
1733
1734 @heading Disadvantages
1735
1736 Of course it also has some disadvantages:
1737
1738 @itemize @bullet
1739
1740 @item
1741 not interactive: GiNaC programs have to be written in 
1742 an editor, compiled and executed. You cannot play with 
1743 expressions interactively.  However, such an extension is not
1744 inherently forbidden by design.  In fact, two interactive
1745 interfaces are possible: First, a simple shell that exposes GiNaC's
1746 types to a command line can readily be written (and has been
1747 written) and second, as a more consistent approach we plan
1748 an integration with the @acronym{CINT} C++ interpreter.
1749
1750 @item
1751 advanced features: GiNaC cannot compete with a program
1752 like @emph{Reduce} which exists for more than
1753 30 years now or @emph{Maple} which grows since 
1754 1981 by the work of dozens of programmers, with respect to
1755 mathematical features. Integration, factorization, non-trivial
1756 simplifications, limits etc. are missing in GiNaC (and are not
1757 planned for the near future).
1758
1759 @item
1760 portability: While the GiNaC library itself is designed
1761 to avoid any platform dependent features (it should compile
1762 on any ANSI compliant C++ compiler), the
1763 currently used version of the CLN library (fast large integer and
1764 arbitrary precision arithmetics) can be compiled only on systems
1765 with a recently new C++ compiler from the
1766 GNU Compiler Collection (@acronym{GCC}).  GiNaC uses
1767 recent language features like explicit constructors, mutable
1768 members, RTTI, @code{dynamic_cast}s and STL, so ANSI compliance is meant
1769 literally.  Recent @acronym{GCC} versions starting at
1770 2.95, although itself not yet ANSI compliant, support all needed
1771 features.
1772     
1773 @end itemize
1774
1775
1776 @heading Why C++?
1777
1778 Why did we choose to implement GiNaC in C++ instead of Java or any other
1779 language? C++ is not perfect: type checking is not strict
1780 (casting is possible), separation between interface and implementation
1781 is not complete, object oriented design is not enforced.  The main
1782 reason is the often scolded feature of operator overloading in
1783 C++. While it may be true that operating on classes
1784 with a @code{+} operator is rarely meaningful, it is
1785 perfectly suited for algebraic expressions. Writing @math{3x+5y} as
1786 @code{3*x+5*y} instead of @code{x.times(3).plus(y.times(5))} looks much more
1787 natural. Furthermore, the main developers are more familiar with
1788 C++ than with any other programming language.
1789
1790
1791 @node Bibliography, Concept Index, A Comparison With Other CAS, Top
1792 @c    node-name, next, previous, up
1793 @chapter Bibliography
1794
1795 @itemize @minus{}
1796
1797 @item
1798 @cite{ISO/IEC 14882:1998: Programming Languages: C++}
1799
1800 @item
1801 @cite{CLN: A Class Library for Numbers}, @email{haible@@ilog.fr, Bruno Haible}
1802
1803 @item
1804 @cite{The C++ Programming Language}, Bjarne Stroustrup, 3rd Edition, ISBN 0-201-88954-4, Addison Wesley
1805
1806 @item
1807 @cite{C++ FAQs}, Marshall Cline, ISBN 0-201-58958-3, 1995, Addison Wesley
1808
1809 @item
1810 @cite{Algorithms for Computer Algebra}, Keith O. Geddes, Stephen R. Czapor,
1811 and George Labahn, ISBN 0-7923-9259-0, 1992, Kluwer Academic Publishers, Norwell, Massachusetts
1812
1813 @item
1814 @cite{Computer Algebra: Systems and Algorithms for Algebraic Computation},
1815 J.H. Davenport, Y. Siret, and E. Tournier, ISBN 0-12-204230-1, 1988, 
1816 Academic Press, London
1817
1818 @end itemize
1819
1820
1821 @node Concept Index, , Bibliography, Top
1822 @c    node-name, next, previous, up
1823 @unnumbered Concept Index
1824
1825 @printindex cp
1826
1827 @bye
1828