ginac/basic.cpp

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