ginac/ex.cpp

   1 /** @file ex.cpp
   2  *
   3  *  Implementation of GiNaC's light-weight expression handles. */
   4
   5 #include <iostream>
   6
   7 #include "ginac.h"
   8
   9 //////////
  10 // default constructor, destructor, copy constructor assignment operator and helpers
  11 //////////
  12
  13 // public
  14
  15 #ifndef INLINE_EX_CONSTRUCTORS
  16
  17 ex::ex() : bp(exZERO().bp)
  18 {
  19     debugmsg("ex default constructor",LOGLEVEL_CONSTRUCT);
  20     ASSERT(exZERO().bp!=0);
  21     ASSERT(exZERO().bp->flags & status_flags::dynallocated);
  22     ASSERT(bp!=0);
  23     ++bp->refcount;
  24 }
  25
  26 ex::~ex()
  27 {
  28     debugmsg("ex destructor",LOGLEVEL_DESTRUCT);
  29     ASSERT(bp!=0);
  30     ASSERT(bp->flags & status_flags::dynallocated);
  31     if (--bp->refcount == 0) {
  32         delete bp;
  33     }
  34 }
  35
  36 ex::ex(ex const & other) : bp(other.bp)
  37 {
  38     debugmsg("ex copy constructor",LOGLEVEL_CONSTRUCT);
  39     ASSERT(bp!=0);
  40     ASSERT((bp->flags) & status_flags::dynallocated);
  41     ++bp->refcount;
  42 }
  43
  44 ex const & ex::operator=(ex const & other)
  45 {
  46     debugmsg("ex operator=",LOGLEVEL_ASSIGNMENT);
  47     ASSERT(bp!=0);
  48     ASSERT(bp->flags & status_flags::dynallocated);
  49     ASSERT(other.bp!=0);
  50     ASSERT(other.bp->flags & status_flags::dynallocated);
  51     ++other.bp->refcount;
  52     basic * tmpbp=other.bp;
  53     if (--bp->refcount==0) {
  54             delete bp;
  55     }
  56     bp=tmpbp;
  57     return *this;
  58 }
  59
  60 #endif // ndef INLINE_EX_CONSTRUCTORS
  61
  62 //////////
  63 // other constructors
  64 //////////
  65
  66 // public
  67
  68 #ifndef INLINE_EX_CONSTRUCTORS
  69 ex::ex(basic const & other)
  70 {
  71     debugmsg("ex constructor from basic",LOGLEVEL_CONSTRUCT);
  72     construct_from_basic(other);
  73 }
  74 #endif
  75
  76 ex::ex(int const i)
  77 {
  78     debugmsg("ex constructor from int",LOGLEVEL_CONSTRUCT);
  79     construct_from_basic(numeric(i));
  80 }
  81
  82 ex::ex(unsigned int const i)
  83 {
  84     debugmsg("ex constructor from unsigned int",LOGLEVEL_CONSTRUCT);
  85     construct_from_basic(numeric(i));
  86 }
  87
  88 ex::ex(long const i)
  89 {
  90     debugmsg("ex constructor from long",LOGLEVEL_CONSTRUCT);
  91     construct_from_basic(numeric(i));
  92 }
  93
  94 ex::ex(unsigned long const i)
  95 {
  96     debugmsg("ex constructor from unsigned long",LOGLEVEL_CONSTRUCT);
  97     construct_from_basic(numeric(i));
  98 }
  99
 100 ex::ex(double const d)
 101 {
 102     debugmsg("ex constructor from double",LOGLEVEL_CONSTRUCT);
 103     construct_from_basic(numeric(d));
 104 }
 105
 106 //////////
 107 // functions overriding virtual functions from bases classes
 108 //////////
 109
 110 // none
 111
 112 //////////
 113 // new virtual functions which can be overridden by derived classes
 114 //////////
 115
 116 // none
 117
 118 //////////
 119 // non-virtual functions in this class
 120 //////////
 121
 122 // public
 123
 124 void ex::swap(ex & other)
 125 {
 126     debugmsg("ex swap",LOGLEVEL_MEMBER_FUNCTION);
 127
 128     ASSERT(bp!=0);
 129     ASSERT(bp->flags & status_flags::dynallocated);
 130     ASSERT(other.bp!=0);
 131     ASSERT(other.bp->flags & status_flags::dynallocated);
 132
 133     basic * tmpbp=bp;
 134     bp=other.bp;
 135     other.bp=tmpbp;
 136 }
 137
 138 bool ex::info(unsigned inf) const
 139 {
 140     if (inf == info_flags::normal_form) {
 141
 142     // Polynomials are in normal form
 143     if (info(info_flags::polynomial))
 144         return true;
 145
 146     // polynomial^(-int) is in normal form
 147     if (is_ex_exactly_of_type(*this, power))
 148         return op(1).info(info_flags::negint);
 149
 150     // polynomial^(int) * polynomial^(int) * ... is in normal form
 151     if (!is_ex_exactly_of_type(*this, mul))
 152         return false;
 153     for (int i=0; i<nops(); i++) {
 154         if (is_ex_exactly_of_type(op(i), power)) {
 155             if (!op(i).op(1).info(info_flags::integer))
 156                 return false;
 157             if (!op(i).op(0).info(info_flags::polynomial))
 158                 return false;
 159         } else
 160             if (!op(i).info(info_flags::polynomial))
 161                 return false;
 162     }
 163     return true;
 164     } else {
 165         return bp->info(inf);
 166     }
 167 }
 168
 169 int ex::nops() const
 170 {
 171     ASSERT(bp!=0);
 172     return bp->nops();
 173 }
 174
 175 ex ex::expand(unsigned options) const
 176 {
 177     ASSERT(bp!=0);
 178     return bp->expand(options);
 179 }
 180
 181 bool ex::has(ex const & other) const
 182 {
 183     ASSERT(bp!=0);
 184     return bp->has(other);
 185 }
 186
 187 int ex::degree(symbol const & s) const
 188 {
 189     ASSERT(bp!=0);
 190     return bp->degree(s);
 191 }
 192
 193 int ex::ldegree(symbol const & s) const
 194 {
 195     ASSERT(bp!=0);
 196     return bp->ldegree(s);
 197 }
 198
 199 ex ex::coeff(symbol const & s, int const n) const
 200 {
 201     ASSERT(bp!=0);
 202     return bp->coeff(s,n);
 203 }
 204
 205 ex ex::numer(bool normalize) const
 206 {
 207     ex n;
 208     if (normalize && !info(info_flags::normal_form))
 209         n = normal();
 210     else
 211         n = *this;
 212
 213     // polynomial
 214     if (n.info(info_flags::polynomial))
 215         return n;
 216
 217     // something^(-int)
 218     if (is_ex_exactly_of_type(n, power) && n.op(1).info(info_flags::negint))
 219         return exONE();
 220
 221     // something^(int) * something^(int) * ...
 222     if (!is_ex_exactly_of_type(n, mul))
 223         return n;
 224     ex res = exONE();
 225     for (int i=0; i<n.nops(); i++) {
 226         if (!is_ex_exactly_of_type(n.op(i), power) || !n.op(i).op(1).info(info_flags::negint))
 227             res *= n.op(i);
 228     }
 229     return res;
 230 }
 231
 232 ex ex::denom(bool normalize) const
 233 {
 234     ex n;
 235     if (normalize && !info(info_flags::normal_form))
 236         n = normal();
 237     else
 238         n = *this;
 239
 240     // polynomial
 241     if (n.info(info_flags::polynomial))
 242         return exONE();
 243
 244     // something^(-int)
 245     if (is_ex_exactly_of_type(n, power) && n.op(1).info(info_flags::negint))
 246         return power(n.op(0), -(n.op(1)));
 247
 248     // something^(int) * something^(int) * ...
 249     if (!is_ex_exactly_of_type(n, mul))
 250         return exONE();
 251     ex res = exONE();
 252     for (int i=0; i<n.nops(); i++) {
 253         if (is_ex_exactly_of_type(n.op(i), power) && n.op(i).op(1).info(info_flags::negint))
 254             res *= power(n.op(i), -1);
 255     }
 256     return res;
 257 }
 258
 259 ex ex::collect(symbol const & s) const
 260 {
 261     ASSERT(bp!=0);
 262     return bp->collect(s);
 263 }
 264
 265 ex ex::eval(int level) const
 266 {
 267     ASSERT(bp!=0);
 268     return bp->eval(level);
 269 }
 270
 271 ex ex::evalf(int level) const
 272 {
 273     ASSERT(bp!=0);
 274     return bp->evalf(level);
 275 }
 276
 277 ex ex::subs(lst const & ls, lst const & lr) const
 278 {
 279     ASSERT(bp!=0);
 280     return bp->subs(ls,lr);
 281 }
 282
 283 ex ex::subs(ex const & e) const
 284 {
 285     ASSERT(bp!=0);
 286     return bp->subs(e);
 287 }
 288
 289 exvector ex::get_indices(void) const
 290 {
 291     ASSERT(bp!=0);
 292     return bp->get_indices();
 293 }
 294
 295 ex ex::simplify_ncmul(exvector const & v) const
 296 {
 297     ASSERT(bp!=0);
 298     return bp->simplify_ncmul(v);
 299 }
 300
 301 ex ex::operator[](ex const & index) const
 302 {
 303     debugmsg("ex operator[ex]",LOGLEVEL_OPERATOR);
 304     ASSERT(bp!=0);
 305     return (*bp)[index];
 306 }
 307
 308 ex ex::operator[](int const i) const
 309 {
 310     debugmsg("ex operator[int]",LOGLEVEL_OPERATOR);
 311     ASSERT(bp!=0);
 312     return (*bp)[i];
 313 }
 314
 315 ex ex::op(int const i) const
 316 {
 317     debugmsg("ex op()",LOGLEVEL_MEMBER_FUNCTION);
 318     ASSERT(bp!=0);
 319     return bp->op(i);
 320 }
 321
 322 ex & ex::let_op(int const i)
 323 {
 324     debugmsg("ex let_op()",LOGLEVEL_MEMBER_FUNCTION);
 325     makewriteable();
 326     ASSERT(bp!=0);
 327     return bp->let_op(i);
 328 }
 329
 330 #ifndef INLINE_EX_CONSTRUCTORS
 331 int ex::compare(ex const & other) const
 332 {
 333     ASSERT(bp!=0);
 334     ASSERT(other.bp!=0);
 335     if (bp==other.bp) {
 336         // special case: both expression point to same basic, trivially equal
 337         return 0;
 338     }
 339     return bp->compare(*other.bp);
 340 }
 341 #endif // ndef INLINE_EX_CONSTRUCTORS
 342
 343 #ifndef INLINE_EX_CONSTRUCTORS
 344 bool ex::is_equal(ex const & other) const
 345 {
 346     ASSERT(bp!=0);
 347     ASSERT(other.bp!=0);
 348     if (bp==other.bp) {
 349         // special case: both expression point to same basic, trivially equal
 350         return true;
 351     }
 352     return bp->is_equal(*other.bp);
 353 }
 354 #endif // ndef INLINE_EX_CONSTRUCTORS
 355
 356 unsigned ex::return_type(void) const
 357 {
 358     ASSERT(bp!=0);
 359     return bp->return_type();
 360 }
 361
 362 unsigned ex::return_type_tinfo(void) const
 363 {
 364     ASSERT(bp!=0);
 365     return bp->return_type_tinfo();
 366 }
 367
 368 unsigned ex::gethash(void) const
 369 {
 370     ASSERT(bp!=0);
 371     return bp->gethash();
 372 }
 373
 374 ex ex::exadd(ex const & rh) const
 375 {
 376     return (new add(*this,rh))->setflag(status_flags::dynallocated);
 377 }
 378
 379 ex ex::exmul(ex const & rh) const
 380 {
 381     return (new mul(*this,rh))->setflag(status_flags::dynallocated);
 382 }
 383
 384 ex ex::exncmul(ex const & rh) const
 385 {
 386     return (new ncmul(*this,rh))->setflag(status_flags::dynallocated);
 387 }
 388
 389 // private
 390
 391 void ex::makewriteable()
 392 {
 393     debugmsg("ex makewriteable",LOGLEVEL_MEMBER_FUNCTION);
 394     ASSERT(bp!=0);
 395     ASSERT(bp->flags & status_flags::dynallocated);
 396     if (bp->refcount > 1) {
 397         basic * bp2=bp->duplicate();
 398         ++bp2->refcount;
 399         bp2->setflag(status_flags::dynallocated);
 400         --bp->refcount;
 401         bp=bp2;
 402     }
 403     ASSERT(bp->refcount == 1);
 404 }
 405
 406 void ex::construct_from_basic(basic const & other)
 407 {
 408     if ( (other.flags & status_flags::evaluated)==0 ) {
 409         // cf. copy constructor
 410         ex const & tmpex = other.eval(1); // evaluate only one (top) level
 411         bp = tmpex.bp;
 412         ASSERT(bp!=0);
 413         ASSERT(bp->flags & status_flags::dynallocated);
 414         ++bp->refcount;
 415         if ((other.flags & status_flags::dynallocated)&&(other.refcount==0)) {
 416             delete &const_cast<basic &>(other);
 417         }
 418     } else {
 419         if (other.flags & status_flags::dynallocated) {
 420             bp=&const_cast<basic &>(other);
 421         } else {
 422             bp=other.duplicate();
 423             bp->setflag(status_flags::dynallocated);
 424         }
 425         ASSERT(bp!=0);
 426         // bp->clearflag(status_flags::evaluated);
 427         ++bp->refcount;
 428     }
 429     ASSERT(bp!=0);
 430     ASSERT(bp->flags & status_flags::dynallocated);
 431 }
 432
 433 //////////
 434 // static member variables
 435 //////////
 436
 437 // none
 438
 439 //////////
 440 // functions which are not member functions
 441 //////////
 442
 443 // none
 444
 445 //////////
 446 // global functions
 447 //////////
 448
 449 ex const & exZERO(void)
 450 {
 451     static ex * eZERO=new ex(numZERO());
 452     return *eZERO;
 453 }
 454
 455 ex const & exONE(void)
 456 {
 457     static ex * eONE=new ex(numONE());
 458     return *eONE;
 459 }
 460
 461 ex const & exTWO(void)
 462 {
 463     static ex * eTWO=new ex(numTWO());
 464     return *eTWO;
 465 }
 466
 467 ex const & exTHREE(void)
 468 {
 469     static ex * eTHREE=new ex(numTHREE());
 470     return *eTHREE;
 471 }
 472
 473 ex const & exMINUSONE(void)
 474 {
 475     static ex * eMINUSONE=new ex(numMINUSONE());
 476     return *eMINUSONE;
 477 }
 478
 479 ex const & exHALF(void)
 480 {
 481     static ex * eHALF=new ex(ex(1)/ex(2));
 482     return *eHALF;
 483 }
 484
 485 ex const & exMINUSHALF(void)
 486 {
 487     static ex * eMINUSHALF=new ex(numeric(-1,2));
 488     return *eMINUSHALF;
 489 }
 490