implemented double dispatch for print(); methods are specified with
[ginac.git] / ginac / basic.cpp
1 /** @file basic.cpp
2  *
3  *  Implementation of GiNaC's ABC. */
4
5 /*
6  *  GiNaC Copyright (C) 1999-2003 Johannes Gutenberg University Mainz, Germany
7  *
8  *  This program is free software; you can redistribute it and/or modify
9  *  it under the terms of the GNU General Public License as published by
10  *  the Free Software Foundation; either version 2 of the License, or
11  *  (at your option) any later version.
12  *
13  *  This program is distributed in the hope that it will be useful,
14  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16  *  GNU General Public License for more details.
17  *
18  *  You should have received a copy of the GNU General Public License
19  *  along with this program; if not, write to the Free Software
20  *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
21  */
22
23 #include <iostream>
24 #include <stdexcept>
25 #ifdef DO_GINAC_ASSERT
26 #  include <typeinfo>
27 #endif
28
29 #include "basic.h"
30 #include "ex.h"
31 #include "numeric.h"
32 #include "power.h"
33 #include "symbol.h"
34 #include "lst.h"
35 #include "ncmul.h"
36 #include "relational.h"
37 #include "operators.h"
38 #include "wildcard.h"
39 #include "archive.h"
40 #include "utils.h"
41
42 namespace GiNaC {
43
44 GINAC_IMPLEMENT_REGISTERED_CLASS_OPT(basic, void,
45   print_func<print_context>(&basic::do_print).
46   print_func<print_tree>(&basic::do_print_tree).
47   print_func<print_python_repr>(&basic::do_print_python_repr))
48
49 //////////
50 // default constructor, destructor, copy constructor and assignment operator
51 //////////
52
53 // public
54
55 /** basic copy constructor: implicitly assumes that the other class is of
56  *  the exact same type (as it's used by duplicate()), so it can copy the
57  *  tinfo_key and the hash value. */
58 basic::basic(const basic & other) : tinfo_key(other.tinfo_key), flags(other.flags & ~status_flags::dynallocated), hashvalue(other.hashvalue), refcount(0)
59 {
60         GINAC_ASSERT(typeid(*this) == typeid(other));
61 }
62
63 /** basic assignment operator: the other object might be of a derived class. */
64 const basic & basic::operator=(const basic & other)
65 {
66         unsigned fl = other.flags & ~status_flags::dynallocated;
67         if (tinfo_key != other.tinfo_key) {
68                 // The other object is of a derived class, so clear the flags as they
69                 // might no longer apply (especially hash_calculated). Oh, and don't
70                 // copy the tinfo_key: it is already set correctly for this object.
71                 flags = 0;
72         } else {
73                 // The objects are of the exact same class, so copy the hash value.
74                 hashvalue = other.hashvalue;
75         }
76         flags = fl;
77         refcount = 0;
78         return *this;
79 }
80
81 // protected
82
83 // none (all inlined)
84
85 //////////
86 // other constructors
87 //////////
88
89 // none (all inlined)
90
91 //////////
92 // archiving
93 //////////
94
95 /** Construct object from archive_node. */
96 basic::basic(const archive_node &n, lst &sym_lst) : flags(0), refcount(0)
97 {
98         // Reconstruct tinfo_key from class name
99         std::string class_name;
100         if (n.find_string("class", class_name))
101                 tinfo_key = find_tinfo_key(class_name);
102         else
103                 throw (std::runtime_error("archive node contains no class name"));
104 }
105
106 /** Unarchive the object. */
107 DEFAULT_UNARCHIVE(basic)
108
109 /** Archive the object. */
110 void basic::archive(archive_node &n) const
111 {
112         n.add_string("class", class_name());
113 }
114
115 //////////
116 // new virtual functions which can be overridden by derived classes
117 //////////
118
119 // public
120
121 /** Output to stream. This performs double dispatch on the dynamic type of
122  *  *this and the dynamic type of the supplied print context.
123  *  @param c print context object that describes the output formatting
124  *  @param level value that is used to identify the precedence or indentation
125  *               level for placing parentheses and formatting */
126 void basic::print(const print_context & c, unsigned level) const
127 {
128         // Double dispatch on object type and print_context type
129         const registered_class_info * reg_info = &get_class_info();
130         const print_context_class_info * pc_info = &c.get_class_info();
131
132 next_class:
133         const std::vector<print_functor> & pdt = reg_info->options.get_print_dispatch_table();
134
135 next_context:
136         unsigned id = pc_info->options.get_id();
137         if (id >= pdt.size() || !(pdt[id].is_valid())) {
138
139                 // Method not found, try parent print_context class
140                 const print_context_class_info * parent_pc_info = pc_info->get_parent();
141                 if (parent_pc_info) {
142                         pc_info = parent_pc_info;
143                         goto next_context;
144                 }
145
146                 // Method still not found, try parent class
147                 const registered_class_info * parent_reg_info = reg_info->get_parent();
148                 if (parent_reg_info) {
149                         reg_info = parent_reg_info;
150                         pc_info = &c.get_class_info();
151                         goto next_class;
152                 }
153
154                 // Method still not found. This shouldn't happen because basic (the
155                 // base class of the algebraic hierarchy) registers a method for
156                 // print_context (the base class of the print context hierarchy),
157                 // so if we end up here, there's something wrong with the class
158                 // registry.
159                 throw (std::runtime_error(std::string("basic::print(): method for ") + class_name() + "/" + c.class_name() + " not found"));
160
161         } else {
162
163                 // Call method
164                 pdt[id](*this, c, level);
165         }
166 }
167
168 /** Default output to stream. */
169 void basic::do_print(const print_context & c, unsigned level) const
170 {
171         c.s << "[" << class_name() << " object]";
172 }
173
174 /** Tree output to stream. */
175 void basic::do_print_tree(const print_tree & c, unsigned level) const
176 {
177         c.s << std::string(level, ' ') << class_name()
178             << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
179             << ", nops=" << nops()
180             << std::endl;
181         for (size_t i=0; i<nops(); ++i)
182                 op(i).print(c, level + c.delta_indent);
183 }
184
185 /** Python parsable output to stream. */
186 void basic::do_print_python_repr(const print_python_repr & c, unsigned level) const
187 {
188         c.s << class_name() << "()";
189 }
190
191 /** Little wrapper around print to be called within a debugger.
192  *  This is needed because you cannot call foo.print(cout) from within the
193  *  debugger because it might not know what cout is.  This method can be
194  *  invoked with no argument and it will simply print to stdout.
195  *
196  *  @see basic::print
197  *  @see basic::dbgprinttree */
198 void basic::dbgprint() const
199 {
200         this->print(std::cerr);
201         std::cerr << std::endl;
202 }
203
204 /** Little wrapper around printtree to be called within a debugger.
205  *
206  *  @see basic::dbgprint */
207 void basic::dbgprinttree() const
208 {
209         this->print(print_tree(std::cerr));
210 }
211
212 /** Return relative operator precedence (for parenthezing output). */
213 unsigned basic::precedence() const
214 {
215         return 70;
216 }
217
218 /** Information about the object.
219  *
220  *  @see class info_flags */
221 bool basic::info(unsigned inf) const
222 {
223         // all possible properties are false for basic objects
224         return false;
225 }
226
227 /** Number of operands/members. */
228 size_t basic::nops() const
229 {
230         // iterating from 0 to nops() on atomic objects should be an empty loop,
231         // and accessing their elements is a range error.  Container objects should
232         // override this.
233         return 0;
234 }
235
236 /** Return operand/member at position i. */
237 ex basic::op(size_t i) const
238 {
239         throw(std::range_error(std::string("basic::op(): ") + class_name() + std::string(" has no operands")));
240 }
241
242 /** Return modifyable operand/member at position i. */
243 ex & basic::let_op(size_t i)
244 {
245         ensure_if_modifiable();
246         throw(std::range_error(std::string("basic::let_op(): ") + class_name() + std::string(" has no operands")));
247 }
248
249 ex basic::operator[](const ex & index) const
250 {
251         if (is_exactly_a<numeric>(index))
252                 return op(static_cast<size_t>(ex_to<numeric>(index).to_int()));
253
254         throw(std::invalid_argument(std::string("non-numeric indices not supported by ") + class_name()));
255 }
256
257 ex basic::operator[](size_t i) const
258 {
259         return op(i);
260 }
261
262 ex & basic::operator[](const ex & index)
263 {
264         if (is_exactly_a<numeric>(index))
265                 return let_op(ex_to<numeric>(index).to_int());
266
267         throw(std::invalid_argument(std::string("non-numeric indices not supported by ") + class_name()));
268 }
269
270 ex & basic::operator[](size_t i)
271 {
272         return let_op(i);
273 }
274
275 /** Test for occurrence of a pattern.  An object 'has' a pattern if it matches
276  *  the pattern itself or one of the children 'has' it.  As a consequence
277  *  (according to the definition of children) given e=x+y+z, e.has(x) is true
278  *  but e.has(x+y) is false. */
279 bool basic::has(const ex & pattern) const
280 {
281         lst repl_lst;
282         if (match(pattern, repl_lst))
283                 return true;
284         for (size_t i=0; i<nops(); i++)
285                 if (op(i).has(pattern))
286                         return true;
287         
288         return false;
289 }
290
291 /** Construct new expression by applying the specified function to all
292  *  sub-expressions (one level only, not recursively). */
293 ex basic::map(map_function & f) const
294 {
295         size_t num = nops();
296         if (num == 0)
297                 return *this;
298
299         basic *copy = duplicate();
300         copy->setflag(status_flags::dynallocated);
301         copy->clearflag(status_flags::hash_calculated | status_flags::expanded);
302         for (size_t i=0; i<num; i++)
303                 copy->let_op(i) = f(copy->op(i));
304         return *copy;
305 }
306
307 /** Return degree of highest power in object s. */
308 int basic::degree(const ex & s) const
309 {
310         return is_equal(ex_to<basic>(s)) ? 1 : 0;
311 }
312
313 /** Return degree of lowest power in object s. */
314 int basic::ldegree(const ex & s) const
315 {
316         return is_equal(ex_to<basic>(s)) ? 1 : 0;
317 }
318
319 /** Return coefficient of degree n in object s. */
320 ex basic::coeff(const ex & s, int n) const
321 {
322         if (is_equal(ex_to<basic>(s)))
323                 return n==1 ? _ex1 : _ex0;
324         else
325                 return n==0 ? *this : _ex0;
326 }
327
328 /** Sort expanded expression in terms of powers of some object(s).
329  *  @param s object(s) to sort in
330  *  @param distributed recursive or distributed form (only used when s is a list) */
331 ex basic::collect(const ex & s, bool distributed) const
332 {
333         ex x;
334         if (is_a<lst>(s)) {
335
336                 // List of objects specified
337                 if (s.nops() == 0)
338                         return *this;
339                 if (s.nops() == 1)
340                         return collect(s.op(0));
341
342                 else if (distributed) {
343
344                         // Get lower/upper degree of all symbols in list
345                         size_t num = s.nops();
346                         struct sym_info {
347                                 ex sym;
348                                 int ldeg, deg;
349                                 int cnt;  // current degree, 'counter'
350                                 ex coeff; // coefficient for degree 'cnt'
351                         };
352                         sym_info *si = new sym_info[num];
353                         ex c = *this;
354                         for (size_t i=0; i<num; i++) {
355                                 si[i].sym = s.op(i);
356                                 si[i].ldeg = si[i].cnt = this->ldegree(si[i].sym);
357                                 si[i].deg = this->degree(si[i].sym);
358                                 c = si[i].coeff = c.coeff(si[i].sym, si[i].cnt);
359                         }
360
361                         while (true) {
362
363                                 // Calculate coeff*x1^c1*...*xn^cn
364                                 ex y = _ex1;
365                                 for (size_t i=0; i<num; i++) {
366                                         int cnt = si[i].cnt;
367                                         y *= power(si[i].sym, cnt);
368                                 }
369                                 x += y * si[num - 1].coeff;
370
371                                 // Increment counters
372                                 size_t n = num - 1;
373                                 while (true) {
374                                         ++si[n].cnt;
375                                         if (si[n].cnt <= si[n].deg) {
376                                                 // Update coefficients
377                                                 ex c;
378                                                 if (n == 0)
379                                                         c = *this;
380                                                 else
381                                                         c = si[n - 1].coeff;
382                                                 for (size_t i=n; i<num; i++)
383                                                         c = si[i].coeff = c.coeff(si[i].sym, si[i].cnt);
384                                                 break;
385                                         }
386                                         if (n == 0)
387                                                 goto done;
388                                         si[n].cnt = si[n].ldeg;
389                                         n--;
390                                 }
391                         }
392
393 done:           delete[] si;
394
395                 } else {
396
397                         // Recursive form
398                         x = *this;
399                         size_t n = s.nops() - 1;
400                         while (true) {
401                                 x = x.collect(s[n]);
402                                 if (n == 0)
403                                         break;
404                                 n--;
405                         }
406                 }
407
408         } else {
409
410                 // Only one object specified
411                 for (int n=this->ldegree(s); n<=this->degree(s); ++n)
412                         x += this->coeff(s,n)*power(s,n);
413         }
414         
415         // correct for lost fractional arguments and return
416         return x + (*this - x).expand();
417 }
418
419 /** Perform automatic non-interruptive term rewriting rules. */
420 ex basic::eval(int level) const
421 {
422         // There is nothing to do for basic objects:
423         return hold();
424 }
425
426 /** Function object to be applied by basic::evalf(). */
427 struct evalf_map_function : public map_function {
428         int level;
429         evalf_map_function(int l) : level(l) {}
430         ex operator()(const ex & e) { return evalf(e, level); }
431 };
432
433 /** Evaluate object numerically. */
434 ex basic::evalf(int level) const
435 {
436         if (nops() == 0)
437                 return *this;
438         else {
439                 if (level == 1)
440                         return *this;
441                 else if (level == -max_recursion_level)
442                         throw(std::runtime_error("max recursion level reached"));
443                 else {
444                         evalf_map_function map_evalf(level - 1);
445                         return map(map_evalf);
446                 }
447         }
448 }
449
450 /** Function object to be applied by basic::evalm(). */
451 struct evalm_map_function : public map_function {
452         ex operator()(const ex & e) { return evalm(e); }
453 } map_evalm;
454
455 /** Evaluate sums, products and integer powers of matrices. */
456 ex basic::evalm() const
457 {
458         if (nops() == 0)
459                 return *this;
460         else
461                 return map(map_evalm);
462 }
463
464 /** Perform automatic symbolic evaluations on indexed expression that
465  *  contains this object as the base expression. */
466 ex basic::eval_indexed(const basic & i) const
467  // this function can't take a "const ex & i" because that would result
468  // in an infinite eval() loop
469 {
470         // There is nothing to do for basic objects
471         return i.hold();
472 }
473
474 /** Add two indexed expressions. They are guaranteed to be of class indexed
475  *  (or a subclass) and their indices are compatible. This function is used
476  *  internally by simplify_indexed().
477  *
478  *  @param self First indexed expression; it's base object is *this
479  *  @param other Second indexed expression
480  *  @return sum of self and other 
481  *  @see ex::simplify_indexed() */
482 ex basic::add_indexed(const ex & self, const ex & other) const
483 {
484         return self + other;
485 }
486
487 /** Multiply an indexed expression with a scalar. This function is used
488  *  internally by simplify_indexed().
489  *
490  *  @param self Indexed expression; it's base object is *this
491  *  @param other Numeric value
492  *  @return product of self and other
493  *  @see ex::simplify_indexed() */
494 ex basic::scalar_mul_indexed(const ex & self, const numeric & other) const
495 {
496         return self * other;
497 }
498
499 /** Try to contract two indexed expressions that appear in the same product. 
500  *  If a contraction exists, the function overwrites one or both of the
501  *  expressions and returns true. Otherwise it returns false. It is
502  *  guaranteed that both expressions are of class indexed (or a subclass)
503  *  and that at least one dummy index has been found. This functions is
504  *  used internally by simplify_indexed().
505  *
506  *  @param self Pointer to first indexed expression; it's base object is *this
507  *  @param other Pointer to second indexed expression
508  *  @param v The complete vector of factors
509  *  @return true if the contraction was successful, false otherwise
510  *  @see ex::simplify_indexed() */
511 bool basic::contract_with(exvector::iterator self, exvector::iterator other, exvector & v) const
512 {
513         // Do nothing
514         return false;
515 }
516
517 /** Check whether the expression matches a given pattern. For every wildcard
518  *  object in the pattern, an expression of the form "wildcard == matching_expression"
519  *  is added to repl_lst. */
520 bool basic::match(const ex & pattern, lst & repl_lst) const
521 {
522 /*
523         Sweet sweet shapes, sweet sweet shapes,
524         That's the key thing, right right.
525         Feed feed face, feed feed shapes,
526         But who is the king tonight?
527         Who is the king tonight?
528         Pattern is the thing, the key thing-a-ling,
529         But who is the king of Pattern?
530         But who is the king, the king thing-a-ling,
531         Who is the king of Pattern?
532         Bog is the king, the king thing-a-ling,
533         Bog is the king of Pattern.
534         Ba bu-bu-bu-bu bu-bu-bu-bu-bu-bu bu-bu
535         Bog is the king of Pattern.
536 */
537
538         if (is_exactly_a<wildcard>(pattern)) {
539
540                 // Wildcard matches anything, but check whether we already have found
541                 // a match for that wildcard first (if so, it the earlier match must
542                 // be the same expression)
543                 for (lst::const_iterator it = repl_lst.begin(); it != repl_lst.end(); ++it) {
544                         if (it->op(0).is_equal(pattern))
545                                 return is_equal(ex_to<basic>(it->op(1)));
546                 }
547                 repl_lst.append(pattern == *this);
548                 return true;
549
550         } else {
551
552                 // Expression must be of the same type as the pattern
553                 if (tinfo() != ex_to<basic>(pattern).tinfo())
554                         return false;
555
556                 // Number of subexpressions must match
557                 if (nops() != pattern.nops())
558                         return false;
559
560                 // No subexpressions? Then just compare the objects (there can't be
561                 // wildcards in the pattern)
562                 if (nops() == 0)
563                         return is_equal_same_type(ex_to<basic>(pattern));
564
565                 // Check whether attributes that are not subexpressions match
566                 if (!match_same_type(ex_to<basic>(pattern)))
567                         return false;
568
569                 // Otherwise the subexpressions must match one-to-one
570                 for (size_t i=0; i<nops(); i++)
571                         if (!op(i).match(pattern.op(i), repl_lst))
572                                 return false;
573
574                 // Looks similar enough, match found
575                 return true;
576         }
577 }
578
579 /** Helper function for subs(). Does not recurse into subexpressions. */
580 ex basic::subs_one_level(const exmap & m, unsigned options) const
581 {
582         exmap::const_iterator it;
583
584         if (options & subs_options::no_pattern) {
585                 it = m.find(*this);
586                 if (it != m.end())
587                         return it->second;
588         } else {
589                 for (it = m.begin(); it != m.end(); ++it) {
590                         lst repl_lst;
591                         if (match(ex_to<basic>(it->first), repl_lst))
592                                 return it->second.subs(repl_lst, options | subs_options::no_pattern); // avoid infinite recursion when re-substituting the wildcards
593                 }
594         }
595
596         return *this;
597 }
598
599 /** Substitute a set of objects by arbitrary expressions. The ex returned
600  *  will already be evaluated. */
601 ex basic::subs(const exmap & m, unsigned options) const
602 {
603         size_t num = nops();
604         if (num) {
605
606                 // Substitute in subexpressions
607                 for (size_t i=0; i<num; i++) {
608                         const ex & orig_op = op(i);
609                         const ex & subsed_op = orig_op.subs(m, options);
610                         if (!are_ex_trivially_equal(orig_op, subsed_op)) {
611
612                                 // Something changed, clone the object
613                                 basic *copy = duplicate();
614                                 copy->setflag(status_flags::dynallocated);
615                                 copy->clearflag(status_flags::hash_calculated | status_flags::expanded);
616
617                                 // Substitute the changed operand
618                                 copy->let_op(i++) = subsed_op;
619
620                                 // Substitute the other operands
621                                 for (; i<num; i++)
622                                         copy->let_op(i) = op(i).subs(m, options);
623
624                                 // Perform substitutions on the new object as a whole
625                                 return copy->subs_one_level(m, options);
626                         }
627                 }
628         }
629
630         // Nothing changed or no subexpressions
631         return subs_one_level(m, options);
632 }
633
634 /** Default interface of nth derivative ex::diff(s, n).  It should be called
635  *  instead of ::derivative(s) for first derivatives and for nth derivatives it
636  *  just recurses down.
637  *
638  *  @param s symbol to differentiate in
639  *  @param nth order of differentiation
640  *  @see ex::diff */
641 ex basic::diff(const symbol & s, unsigned nth) const
642 {
643         // trivial: zeroth derivative
644         if (nth==0)
645                 return ex(*this);
646         
647         // evaluate unevaluated *this before differentiating
648         if (!(flags & status_flags::evaluated))
649                 return ex(*this).diff(s, nth);
650         
651         ex ndiff = this->derivative(s);
652         while (!ndiff.is_zero() &&    // stop differentiating zeros
653                nth>1) {
654                 ndiff = ndiff.diff(s);
655                 --nth;
656         }
657         return ndiff;
658 }
659
660 /** Return a vector containing the free indices of an expression. */
661 exvector basic::get_free_indices() const
662 {
663         return exvector(); // return an empty exvector
664 }
665
666 ex basic::eval_ncmul(const exvector & v) const
667 {
668         return hold_ncmul(v);
669 }
670
671 // protected
672
673 /** Function object to be applied by basic::derivative(). */
674 struct derivative_map_function : public map_function {
675         const symbol &s;
676         derivative_map_function(const symbol &sym) : s(sym) {}
677         ex operator()(const ex & e) { return diff(e, s); }
678 };
679
680 /** Default implementation of ex::diff(). It maps the operation on the
681  *  operands (or returns 0 when the object has no operands).
682  *
683  *  @see ex::diff */
684 ex basic::derivative(const symbol & s) const
685 {
686         if (nops() == 0)
687                 return _ex0;
688         else {
689                 derivative_map_function map_derivative(s);
690                 return map(map_derivative);
691         }
692 }
693
694 /** Returns order relation between two objects of same type.  This needs to be
695  *  implemented by each class. It may never return anything else than 0,
696  *  signalling equality, or +1 and -1 signalling inequality and determining
697  *  the canonical ordering.  (Perl hackers will wonder why C++ doesn't feature
698  *  the spaceship operator <=> for denoting just this.) */
699 int basic::compare_same_type(const basic & other) const
700 {
701         return compare_pointers(this, &other);
702 }
703
704 /** Returns true if two objects of same type are equal.  Normally needs
705  *  not be reimplemented as long as it wasn't overwritten by some parent
706  *  class, since it just calls compare_same_type().  The reason why this
707  *  function exists is that sometimes it is easier to determine equality
708  *  than an order relation and then it can be overridden. */
709 bool basic::is_equal_same_type(const basic & other) const
710 {
711         return compare_same_type(other)==0;
712 }
713
714 /** Returns true if the attributes of two objects are similar enough for
715  *  a match. This function must not match subexpressions (this is already
716  *  done by basic::match()). Only attributes not accessible by op() should
717  *  be compared. This is also the reason why this function doesn't take the
718  *  wildcard replacement list from match() as an argument: only subexpressions
719  *  are subject to wildcard matches. Also, this function only needs to be
720  *  implemented for container classes because is_equal_same_type() is
721  *  automatically used instead of match_same_type() if nops() == 0.
722  *
723  *  @see basic::match */
724 bool basic::match_same_type(const basic & other) const
725 {
726         // The default is to only consider subexpressions, but not any other
727         // attributes
728         return true;
729 }
730
731 unsigned basic::return_type() const
732 {
733         return return_types::commutative;
734 }
735
736 unsigned basic::return_type_tinfo() const
737 {
738         return tinfo();
739 }
740
741 /** Compute the hash value of an object and if it makes sense to store it in
742  *  the objects status_flags, do so.  The method inherited from class basic
743  *  computes a hash value based on the type and hash values of possible
744  *  members.  For this reason it is well suited for container classes but
745  *  atomic classes should override this implementation because otherwise they
746  *  would all end up with the same hashvalue. */
747 unsigned basic::calchash() const
748 {
749         unsigned v = golden_ratio_hash(tinfo());
750         for (size_t i=0; i<nops(); i++) {
751                 v = rotate_left(v);
752                 v ^= this->op(i).gethash();
753         }
754
755         // store calculated hash value only if object is already evaluated
756         if (flags & status_flags::evaluated) {
757                 setflag(status_flags::hash_calculated);
758                 hashvalue = v;
759         }
760
761         return v;
762 }
763
764 /** Function object to be applied by basic::expand(). */
765 struct expand_map_function : public map_function {
766         unsigned options;
767         expand_map_function(unsigned o) : options(o) {}
768         ex operator()(const ex & e) { return expand(e, options); }
769 };
770
771 /** Expand expression, i.e. multiply it out and return the result as a new
772  *  expression. */
773 ex basic::expand(unsigned options) const
774 {
775         if (nops() == 0)
776                 return (options == 0) ? setflag(status_flags::expanded) : *this;
777         else {
778                 expand_map_function map_expand(options);
779                 return ex_to<basic>(map(map_expand)).setflag(options == 0 ? status_flags::expanded : 0);
780         }
781 }
782
783
784 //////////
785 // non-virtual functions in this class
786 //////////
787
788 // public
789
790 /** Compare objects syntactically to establish canonical ordering.
791  *  All compare functions return: -1 for *this less than other, 0 equal,
792  *  1 greater. */
793 int basic::compare(const basic & other) const
794 {
795         const unsigned hash_this = gethash();
796         const unsigned hash_other = other.gethash();
797         if (hash_this<hash_other) return -1;
798         if (hash_this>hash_other) return 1;
799
800         const unsigned typeid_this = tinfo();
801         const unsigned typeid_other = other.tinfo();
802         if (typeid_this==typeid_other) {
803                 GINAC_ASSERT(typeid(*this)==typeid(other));
804 //              int cmpval = compare_same_type(other);
805 //              if (cmpval!=0) {
806 //                      std::cout << "hash collision, same type: " 
807 //                                << *this << " and " << other << std::endl;
808 //                      this->print(print_tree(std::cout));
809 //                      std::cout << " and ";
810 //                      other.print(print_tree(std::cout));
811 //                      std::cout << std::endl;
812 //              }
813 //              return cmpval;
814                 return compare_same_type(other);
815         } else {
816 //              std::cout << "hash collision, different types: " 
817 //                        << *this << " and " << other << std::endl;
818 //              this->print(print_tree(std::cout));
819 //              std::cout << " and ";
820 //              other.print(print_tree(std::cout));
821 //              std::cout << std::endl;
822                 return (typeid_this<typeid_other ? -1 : 1);
823         }
824 }
825
826 /** Test for syntactic equality.
827  *  This is only a quick test, meaning objects should be in the same domain.
828  *  You might have to .expand(), .normal() objects first, depending on the
829  *  domain of your computation, to get a more reliable answer.
830  *
831  *  @see is_equal_same_type */
832 bool basic::is_equal(const basic & other) const
833 {
834         if (this->gethash()!=other.gethash())
835                 return false;
836         if (this->tinfo()!=other.tinfo())
837                 return false;
838         
839         GINAC_ASSERT(typeid(*this)==typeid(other));
840         
841         return is_equal_same_type(other);
842 }
843
844 // protected
845
846 /** Stop further evaluation.
847  *
848  *  @see basic::eval */
849 const basic & basic::hold() const
850 {
851         return setflag(status_flags::evaluated);
852 }
853
854 /** Ensure the object may be modified without hurting others, throws if this
855  *  is not the case. */
856 void basic::ensure_if_modifiable() const
857 {
858         if (refcount > 1)
859                 throw(std::runtime_error("cannot modify multiply referenced object"));
860         clearflag(status_flags::hash_calculated | status_flags::evaluated);
861 }
862
863 //////////
864 // global variables
865 //////////
866
867 int max_recursion_level = 1024;
868
869 } // namespace GiNaC