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