ginac/mul.cpp

   1 /** @file mul.cpp
   2  *
   3  *  Implementation of GiNaC's products of expressions. */
   4
   5 /*
   6  *  GiNaC Copyright (C) 1999-2000 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 <vector>
  24 #include <stdexcept>
  25
  26 #include "mul.h"
  27 #include "add.h"
  28 #include "power.h"
  29 #include "debugmsg.h"
  30 #include "utils.h"
  31
  32 #ifndef NO_GINAC_NAMESPACE
  33 namespace GiNaC {
  34 #endif // ndef NO_GINAC_NAMESPACE
  35
  36 //////////
  37 // default constructor, destructor, copy constructor assignment operator and helpers
  38 //////////
  39
  40 // public
  41
  42 mul::mul()
  43 {
  44     debugmsg("mul default constructor",LOGLEVEL_CONSTRUCT);
  45     tinfo_key = TINFO_mul;
  46 }
  47
  48 mul::~mul()
  49 {
  50     debugmsg("mul destructor",LOGLEVEL_DESTRUCT);
  51     destroy(0);
  52 }
  53
  54 mul::mul(mul const & other)
  55 {
  56     debugmsg("mul copy constructor",LOGLEVEL_CONSTRUCT);
  57     copy(other);
  58 }
  59
  60 mul const & mul::operator=(mul const & other)
  61 {
  62     debugmsg("mul operator=",LOGLEVEL_ASSIGNMENT);
  63     if (this != &other) {
  64         destroy(1);
  65         copy(other);
  66     }
  67     return *this;
  68 }
  69
  70 // protected
  71
  72 void mul::copy(mul const & other)
  73 {
  74     expairseq::copy(other);
  75 }
  76
  77 void mul::destroy(bool call_parent)
  78 {
  79     if (call_parent) expairseq::destroy(call_parent);
  80 }
  81
  82 //////////
  83 // other constructors
  84 //////////
  85
  86 // public
  87
  88 mul::mul(ex const & lh, ex const & rh)
  89 {
  90     debugmsg("mul constructor from ex,ex",LOGLEVEL_CONSTRUCT);
  91     tinfo_key = TINFO_mul;
  92     overall_coeff=_ex1();
  93     construct_from_2_ex(lh,rh);
  94     GINAC_ASSERT(is_canonical());
  95 }
  96
  97 mul::mul(exvector const & v)
  98 {
  99     debugmsg("mul constructor from exvector",LOGLEVEL_CONSTRUCT);
 100     tinfo_key = TINFO_mul;
 101     overall_coeff=_ex1();
 102     construct_from_exvector(v);
 103     GINAC_ASSERT(is_canonical());
 104 }
 105
 106 /*
 107 mul::mul(epvector const & v, bool do_not_canonicalize)
 108 {
 109     debugmsg("mul constructor from epvector,bool",LOGLEVEL_CONSTRUCT);
 110     tinfo_key = TINFO_mul;
 111     if (do_not_canonicalize) {
 112         seq=v;
 113 #ifdef EXPAIRSEQ_USE_HASHTAB
 114         combine_same_terms(); // to build hashtab
 115 #endif // def EXPAIRSEQ_USE_HASHTAB
 116     } else {
 117         construct_from_epvector(v);
 118     }
 119     GINAC_ASSERT(is_canonical());
 120 }
 121 */
 122
 123 mul::mul(epvector const & v)
 124 {
 125     debugmsg("mul constructor from epvector",LOGLEVEL_CONSTRUCT);
 126     tinfo_key = TINFO_mul;
 127     overall_coeff=_ex1();
 128     construct_from_epvector(v);
 129     GINAC_ASSERT(is_canonical());
 130 }
 131
 132 mul::mul(epvector const & v, ex const & oc)
 133 {
 134     debugmsg("mul constructor from epvector,ex",LOGLEVEL_CONSTRUCT);
 135     tinfo_key = TINFO_mul;
 136     overall_coeff=oc;
 137     construct_from_epvector(v);
 138     GINAC_ASSERT(is_canonical());
 139 }
 140
 141 mul::mul(epvector * vp, ex const & oc)
 142 {
 143     debugmsg("mul constructor from epvector *,ex",LOGLEVEL_CONSTRUCT);
 144     tinfo_key = TINFO_mul;
 145     GINAC_ASSERT(vp!=0);
 146     overall_coeff=oc;
 147     construct_from_epvector(*vp);
 148     delete vp;
 149     GINAC_ASSERT(is_canonical());
 150 }
 151
 152 mul::mul(ex const & lh, ex const & mh, ex const & rh)
 153 {
 154     debugmsg("mul constructor from ex,ex,ex",LOGLEVEL_CONSTRUCT);
 155     tinfo_key = TINFO_mul;
 156     exvector factors;
 157     factors.reserve(3);
 158     factors.push_back(lh);
 159     factors.push_back(mh);
 160     factors.push_back(rh);
 161     overall_coeff=_ex1();
 162     construct_from_exvector(factors);
 163     GINAC_ASSERT(is_canonical());
 164 }
 165
 166 //////////
 167 // functions overriding virtual functions from bases classes
 168 //////////
 169
 170 // public
 171
 172 basic * mul::duplicate() const
 173 {
 174     debugmsg("mul duplicate",LOGLEVEL_ASSIGNMENT);
 175     return new mul(*this);
 176 }
 177
 178 void mul::print(ostream & os, unsigned upper_precedence) const
 179 {
 180     debugmsg("mul print",LOGLEVEL_PRINT);
 181     if (precedence<=upper_precedence) os << "(";
 182     bool first=true;
 183     // first print the overall numeric coefficient:
 184     numeric coeff = ex_to_numeric(overall_coeff);
 185     if (coeff.csgn()==-1) os << '-';
 186     if (!coeff.is_equal(_num1()) &&
 187         !coeff.is_equal(_num_1())) {
 188         if (coeff.is_rational()) {
 189             if (coeff.is_negative())
 190                 os << -coeff;
 191             else
 192                 os << coeff;
 193         } else {
 194             if (coeff.csgn()==-1)
 195                 (-coeff).print(os, precedence);
 196             else
 197                 coeff.print(os, precedence);
 198         }
 199         os << '*';
 200     }
 201     // then proceed with the remaining factors:
 202     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
 203         if (!first) {
 204             os << '*';
 205         } else {
 206             first=false;
 207         }
 208         recombine_pair_to_ex(*cit).print(os,precedence);
 209     }
 210     if (precedence<=upper_precedence) os << ")";
 211 }
 212
 213 void mul::printraw(ostream & os) const
 214 {
 215     debugmsg("mul printraw",LOGLEVEL_PRINT);
 216
 217     os << "*(";
 218     for (epvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) {
 219         os << "(";
 220         (*it).rest.bp->printraw(os);
 221         os << ",";
 222         (*it).coeff.bp->printraw(os);
 223         os << "),";
 224     }
 225     os << ",hash=" << hashvalue << ",flags=" << flags;
 226     os << ")";
 227 }
 228
 229 void mul::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const
 230 {
 231     debugmsg("mul print csrc", LOGLEVEL_PRINT);
 232     if (precedence <= upper_precedence)
 233         os << "(";
 234
 235     if (!overall_coeff.is_equal(_ex1())) {
 236         overall_coeff.bp->printcsrc(os,type,precedence);
 237         os << "*";
 238     }
 239
 240     // Print arguments, separated by "*" or "/"
 241     epvector::const_iterator it = seq.begin();
 242     epvector::const_iterator itend = seq.end();
 243     while (it != itend) {
 244
 245         // If the first argument is a negative integer power, it gets printed as "1.0/<expr>"
 246         if (it == seq.begin() && ex_to_numeric(it->coeff).is_integer() && it->coeff.compare(_num0()) < 0) {
 247             if (type == csrc_types::ctype_cl_N)
 248                 os << "recip(";
 249             else
 250                 os << "1.0/";
 251         }
 252
 253         // If the exponent is 1 or -1, it is left out
 254         if (it->coeff.compare(_ex1()) == 0 || it->coeff.compare(_num_1()) == 0)
 255             it->rest.bp->printcsrc(os, type, precedence);
 256         else
 257             // outer parens around ex needed for broken gcc-2.95 parser:
 258             (ex(power(it->rest, abs(ex_to_numeric(it->coeff))))).bp->printcsrc(os, type, upper_precedence);
 259
 260         // Separator is "/" for negative integer powers, "*" otherwise
 261         it++;
 262         if (it != itend) {
 263             if (ex_to_numeric(it->coeff).is_integer() && it->coeff.compare(_num0()) < 0)
 264                 os << "/";
 265             else
 266                 os << "*";
 267         }
 268     }
 269     if (precedence <= upper_precedence)
 270         os << ")";
 271 }
 272
 273 bool mul::info(unsigned inf) const
 274 {
 275     // TODO: optimize
 276     if (inf==info_flags::polynomial ||
 277         inf==info_flags::integer_polynomial ||
 278         inf==info_flags::cinteger_polynomial ||
 279         inf==info_flags::rational_polynomial ||
 280         inf==info_flags::crational_polynomial ||
 281         inf==info_flags::rational_function) {
 282         for (epvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) {
 283             if (!(recombine_pair_to_ex(*it).info(inf)))
 284                 return false;
 285         }
 286         return overall_coeff.info(inf);
 287     } else {
 288         return expairseq::info(inf);
 289     }
 290 }
 291
 292 typedef vector<int> intvector;
 293
 294 int mul::degree(symbol const & s) const
 295 {
 296     int deg_sum=0;
 297     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
 298         deg_sum+=(*cit).rest.degree(s) * ex_to_numeric((*cit).coeff).to_int();
 299     }
 300     return deg_sum;
 301 }
 302
 303 int mul::ldegree(symbol const & s) const
 304 {
 305     int deg_sum=0;
 306     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
 307         deg_sum+=(*cit).rest.ldegree(s) * ex_to_numeric((*cit).coeff).to_int();
 308     }
 309     return deg_sum;
 310 }
 311
 312 ex mul::coeff(symbol const & s, int const n) const
 313 {
 314     exvector coeffseq;
 315     coeffseq.reserve(seq.size()+1);
 316
 317     if (n==0) {
 318         // product of individual coeffs
 319         // if a non-zero power of s is found, the resulting product will be 0
 320         epvector::const_iterator it=seq.begin();
 321         while (it!=seq.end()) {
 322             coeffseq.push_back(recombine_pair_to_ex(*it).coeff(s,n));
 323             ++it;
 324         }
 325         coeffseq.push_back(overall_coeff);
 326         return (new mul(coeffseq))->setflag(status_flags::dynallocated);
 327     }
 328
 329     epvector::const_iterator it=seq.begin();
 330     bool coeff_found=0;
 331     while (it!=seq.end()) {
 332         ex t=recombine_pair_to_ex(*it);
 333         ex c=t.coeff(s,n);
 334         if (!c.is_zero()) {
 335             coeffseq.push_back(c);
 336             coeff_found=1;
 337         } else {
 338             coeffseq.push_back(t);
 339         }
 340         ++it;
 341     }
 342     if (coeff_found) {
 343         coeffseq.push_back(overall_coeff);
 344         return (new mul(coeffseq))->setflag(status_flags::dynallocated);
 345     }
 346
 347     return _ex0();
 348 }
 349
 350 ex mul::eval(int level) const
 351 {
 352     // simplifications  *(...,x;0) -> 0
 353     //                  *(+(x,y,...);c) -> *(+(*(x,c),*(y,c),...)) (c numeric())
 354     //                  *(x;1) -> x
 355     //                  *(;c) -> c
 356
 357     debugmsg("mul eval",LOGLEVEL_MEMBER_FUNCTION);
 358
 359     epvector * evaled_seqp=evalchildren(level);
 360     if (evaled_seqp!=0) {
 361         // do more evaluation later
 362         return (new mul(evaled_seqp,overall_coeff))->
 363                    setflag(status_flags::dynallocated);
 364     }
 365
 366 #ifdef DO_GINAC_ASSERT
 367     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
 368         GINAC_ASSERT((!is_ex_exactly_of_type((*cit).rest,mul))||
 369                (!(ex_to_numeric((*cit).coeff).is_integer())));
 370         GINAC_ASSERT(!((*cit).is_numeric_with_coeff_1()));
 371         if (is_ex_exactly_of_type(recombine_pair_to_ex(*cit),numeric)) {
 372             printtree(cerr,0);
 373         }
 374         GINAC_ASSERT(!is_ex_exactly_of_type(recombine_pair_to_ex(*cit),numeric));
 375         /* for paranoia */
 376         expair p=split_ex_to_pair(recombine_pair_to_ex(*cit));
 377         GINAC_ASSERT(p.rest.is_equal((*cit).rest));
 378         GINAC_ASSERT(p.coeff.is_equal((*cit).coeff));
 379         /* end paranoia */
 380     }
 381 #endif // def DO_GINAC_ASSERT
 382
 383     if (flags & status_flags::evaluated) {
 384         GINAC_ASSERT(seq.size()>0);
 385         GINAC_ASSERT((seq.size()>1)||!overall_coeff.is_equal(_ex1()));
 386         return *this;
 387     }
 388
 389     int seq_size=seq.size();
 390     if (overall_coeff.is_equal(_ex0())) {
 391         // *(...,x;0) -> 0
 392         return _ex0();
 393     } else if (seq_size==0) {
 394         // *(;c) -> c
 395         return overall_coeff;
 396     } else if ((seq_size==1)&&overall_coeff.is_equal(_ex1())) {
 397         // *(x;1) -> x
 398         return recombine_pair_to_ex(*(seq.begin()));
 399     } else if ((seq_size==1) &&
 400                is_ex_exactly_of_type((*seq.begin()).rest,add) &&
 401                ex_to_numeric((*seq.begin()).coeff).is_equal(_num1())) {
 402         // *(+(x,y,...);c) -> +(*(x,c),*(y,c),...) (c numeric(), no powers of +())
 403         add const & addref=ex_to_add((*seq.begin()).rest);
 404         epvector distrseq;
 405         distrseq.reserve(addref.seq.size());
 406         for (epvector::const_iterator cit=addref.seq.begin(); cit!=addref.seq.end(); ++cit) {
 407             distrseq.push_back(addref.combine_pair_with_coeff_to_pair(*cit,
 408                                    overall_coeff));
 409         }
 410         return (new add(distrseq,
 411                         ex_to_numeric(addref.overall_coeff).
 412                         mul_dyn(ex_to_numeric(overall_coeff))))
 413             ->setflag(status_flags::dynallocated  |
 414                       status_flags::evaluated );
 415     }
 416     return this->hold();
 417 }
 418
 419 exvector mul::get_indices(void) const
 420 {
 421     // return union of indices of factors
 422     exvector iv;
 423     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
 424         exvector subiv=(*cit).rest.get_indices();
 425         iv.reserve(iv.size()+subiv.size());
 426         for (exvector::const_iterator cit2=subiv.begin(); cit2!=subiv.end(); ++cit2) {
 427             iv.push_back(*cit2);
 428         }
 429     }
 430     return iv;
 431 }
 432
 433 ex mul::simplify_ncmul(exvector const & v) const
 434 {
 435     throw(std::logic_error("mul::simplify_ncmul() should never have been called!"));
 436 }
 437
 438 // protected
 439
 440 int mul::compare_same_type(basic const & other) const
 441 {
 442     return expairseq::compare_same_type(other);
 443 }
 444
 445 bool mul::is_equal_same_type(basic const & other) const
 446 {
 447     return expairseq::is_equal_same_type(other);
 448 }
 449
 450 unsigned mul::return_type(void) const
 451 {
 452     if (seq.size()==0) {
 453         // mul without factors: should not happen, but commutes
 454         return return_types::commutative;
 455     }
 456
 457     bool all_commutative=1;
 458     unsigned rt;
 459     epvector::const_iterator cit_noncommutative_element; // point to first found nc element
 460
 461     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
 462         rt=(*cit).rest.return_type();
 463         if (rt==return_types::noncommutative_composite) return rt; // one ncc -> mul also ncc
 464         if ((rt==return_types::noncommutative)&&(all_commutative)) {
 465             // first nc element found, remember position
 466             cit_noncommutative_element=cit;
 467             all_commutative=0;
 468         }
 469         if ((rt==return_types::noncommutative)&&(!all_commutative)) {
 470             // another nc element found, compare type_infos
 471             if ((*cit_noncommutative_element).rest.return_type_tinfo()!=(*cit).rest.return_type_tinfo()) {
 472                 // diffent types -> mul is ncc
 473                 return return_types::noncommutative_composite;
 474             }
 475         }
 476     }
 477     // all factors checked
 478     return all_commutative ? return_types::commutative : return_types::noncommutative;
 479 }
 480
 481 unsigned mul::return_type_tinfo(void) const
 482 {
 483     if (seq.size()==0) {
 484         // mul without factors: should not happen
 485         return tinfo_key;
 486     }
 487     // return type_info of first noncommutative element
 488     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
 489         if ((*cit).rest.return_type()==return_types::noncommutative) {
 490             return (*cit).rest.return_type_tinfo();
 491         }
 492     }
 493     // no noncommutative element found, should not happen
 494     return tinfo_key;
 495 }
 496
 497 ex mul::thisexpairseq(epvector const & v, ex const & oc) const
 498 {
 499     return (new mul(v,oc))->setflag(status_flags::dynallocated);
 500 }
 501
 502 ex mul::thisexpairseq(epvector * vp, ex const & oc) const
 503 {
 504     return (new mul(vp,oc))->setflag(status_flags::dynallocated);
 505 }
 506
 507 expair mul::split_ex_to_pair(ex const & e) const
 508 {
 509     if (is_ex_exactly_of_type(e,power)) {
 510         power const & powerref=ex_to_power(e);
 511         if (is_ex_exactly_of_type(powerref.exponent,numeric)) {
 512             return expair(powerref.basis,powerref.exponent);
 513         }
 514     }
 515     return expair(e,_ex1());
 516 }
 517
 518 expair mul::combine_ex_with_coeff_to_pair(ex const & e,
 519                                           ex const & c) const
 520 {
 521     // to avoid duplication of power simplification rules,
 522     // we create a temporary power object
 523     // otherwise it would be hard to correctly simplify
 524     // expression like (4^(1/3))^(3/2)
 525     if (are_ex_trivially_equal(c,_ex1())) {
 526         return split_ex_to_pair(e);
 527     }
 528     return split_ex_to_pair(power(e,c));
 529 }
 530
 531 expair mul::combine_pair_with_coeff_to_pair(expair const & p,
 532                                             ex const & c) const
 533 {
 534     // to avoid duplication of power simplification rules,
 535     // we create a temporary power object
 536     // otherwise it would be hard to correctly simplify
 537     // expression like (4^(1/3))^(3/2)
 538     if (are_ex_trivially_equal(c,_ex1())) {
 539         return p;
 540     }
 541     return split_ex_to_pair(power(recombine_pair_to_ex(p),c));
 542 }
 543
 544 ex mul::recombine_pair_to_ex(expair const & p) const
 545 {
 546     // if (p.coeff.compare(_ex1())==0) {
 547     // if (are_ex_trivially_equal(p.coeff,_ex1())) {
 548     if (ex_to_numeric(p.coeff).is_equal(_num1())) {
 549         return p.rest;
 550     } else {
 551         return power(p.rest,p.coeff);
 552     }
 553 }
 554
 555 bool mul::expair_needs_further_processing(epp it)
 556 {
 557     if (is_ex_exactly_of_type((*it).rest,mul) &&
 558         ex_to_numeric((*it).coeff).is_integer()) {
 559         // combined pair is product with integer power -> expand it
 560         *it=split_ex_to_pair(recombine_pair_to_ex(*it));
 561         return true;
 562     }
 563     if (is_ex_exactly_of_type((*it).rest,numeric)) {
 564         expair ep=split_ex_to_pair(recombine_pair_to_ex(*it));
 565         if (!ep.is_equal(*it)) {
 566             // combined pair is a numeric power which can be simplified
 567             *it=ep;
 568             return true;
 569         }
 570         if (ex_to_numeric((*it).coeff).is_equal(_num1())) {
 571             // combined pair has coeff 1 and must be moved to the end
 572             return true;
 573         }
 574     }
 575     return false;
 576 }
 577
 578 ex mul::default_overall_coeff(void) const
 579 {
 580     return _ex1();
 581 }
 582
 583 void mul::combine_overall_coeff(ex const & c)
 584 {
 585     GINAC_ASSERT(is_ex_exactly_of_type(overall_coeff,numeric));
 586     GINAC_ASSERT(is_ex_exactly_of_type(c,numeric));
 587     overall_coeff = ex_to_numeric(overall_coeff).mul_dyn(ex_to_numeric(c));
 588 }
 589
 590 void mul::combine_overall_coeff(ex const & c1, ex const & c2)
 591 {
 592     GINAC_ASSERT(is_ex_exactly_of_type(overall_coeff,numeric));
 593     GINAC_ASSERT(is_ex_exactly_of_type(c1,numeric));
 594     GINAC_ASSERT(is_ex_exactly_of_type(c2,numeric));
 595     overall_coeff = ex_to_numeric(overall_coeff).
 596                         mul_dyn(ex_to_numeric(c1).power(ex_to_numeric(c2)));
 597 }
 598
 599 bool mul::can_make_flat(expair const & p) const
 600 {
 601     GINAC_ASSERT(is_ex_exactly_of_type(p.coeff,numeric));
 602     // this assertion will probably fail somewhere
 603     // it would require a more careful make_flat, obeying the power laws
 604     // probably should return true only if p.coeff is integer
 605     return ex_to_numeric(p.coeff).is_equal(_num1());
 606 }
 607
 608 ex mul::expand(unsigned options) const
 609 {
 610     exvector sub_expanded_seq;
 611     intvector positions_of_adds;
 612     intvector number_of_add_operands;
 613
 614     epvector * expanded_seqp=expandchildren(options);
 615
 616     epvector const & expanded_seq = expanded_seqp==0 ? seq : *expanded_seqp;
 617
 618     positions_of_adds.resize(expanded_seq.size());
 619     number_of_add_operands.resize(expanded_seq.size());
 620
 621     int number_of_adds=0;
 622     int number_of_expanded_terms=1;
 623
 624     unsigned current_position=0;
 625     epvector::const_iterator last=expanded_seq.end();
 626     for (epvector::const_iterator cit=expanded_seq.begin(); cit!=last; ++cit) {
 627         if (is_ex_exactly_of_type((*cit).rest,add)&&
 628             (ex_to_numeric((*cit).coeff).is_equal(_num1()))) {
 629             positions_of_adds[number_of_adds]=current_position;
 630             add const & expanded_addref=ex_to_add((*cit).rest);
 631             unsigned addref_nops=expanded_addref.nops();
 632             number_of_add_operands[number_of_adds]=addref_nops;
 633             number_of_expanded_terms *= addref_nops;
 634             number_of_adds++;
 635         }
 636         current_position++;
 637     }
 638
 639     if (number_of_adds==0) {
 640         if (expanded_seqp==0) {
 641             return this->setflag(status_flags::expanded);
 642         }
 643         return (new mul(expanded_seqp,overall_coeff))->
 644                      setflag(status_flags::dynallocated ||
 645                              status_flags::expanded);
 646     }
 647
 648     exvector distrseq;
 649     distrseq.reserve(number_of_expanded_terms);
 650
 651     intvector k;
 652     k.resize(number_of_adds);
 653
 654     int l;
 655     for (l=0; l<number_of_adds; l++) {
 656         k[l]=0;
 657     }
 658
 659     while (1) {
 660         epvector term;
 661         term=expanded_seq;
 662         for (l=0; l<number_of_adds; l++) {
 663             add const & addref=ex_to_add(expanded_seq[positions_of_adds[l]].rest);
 664             GINAC_ASSERT(term[positions_of_adds[l]].coeff.compare(_ex1())==0);
 665             term[positions_of_adds[l]]=split_ex_to_pair(addref.op(k[l]));
 666         }
 667         /*
 668         cout << "mul::expand() term begin" << endl;
 669         for (epvector::const_iterator cit=term.begin(); cit!=term.end(); ++cit) {
 670             cout << "rest" << endl;
 671             (*cit).rest.printtree(cout);
 672             cout << "coeff" << endl;
 673             (*cit).coeff.printtree(cout);
 674         }
 675         cout << "mul::expand() term end" << endl;
 676         */
 677         distrseq.push_back((new mul(term,overall_coeff))->
 678                                 setflag(status_flags::dynallocated |
 679                                         status_flags::expanded));
 680
 681         // increment k[]
 682         l=number_of_adds-1;
 683         while ((l>=0)&&((++k[l])>=number_of_add_operands[l])) {
 684             k[l]=0;
 685             l--;
 686         }
 687         if (l<0) break;
 688     }
 689
 690     if (expanded_seqp!=0) {
 691         delete expanded_seqp;
 692     }
 693     /*
 694     cout << "mul::expand() distrseq begin" << endl;
 695     for (exvector::const_iterator cit=distrseq.begin(); cit!=distrseq.end(); ++cit) {
 696         (*cit).printtree(cout);
 697     }
 698     cout << "mul::expand() distrseq end" << endl;
 699     */
 700
 701     return (new add(distrseq))->setflag(status_flags::dynallocated |
 702                                         status_flags::expanded);
 703 }
 704
 705 //////////
 706 // new virtual functions which can be overridden by derived classes
 707 //////////
 708
 709 // none
 710
 711 //////////
 712 // non-virtual functions in this class
 713 //////////
 714
 715 epvector * mul::expandchildren(unsigned options) const
 716 {
 717     epvector::const_iterator last=seq.end();
 718     epvector::const_iterator cit=seq.begin();
 719     while (cit!=last) {
 720         ex const & factor=recombine_pair_to_ex(*cit);
 721         ex const & expanded_factor=factor.expand(options);
 722         if (!are_ex_trivially_equal(factor,expanded_factor)) {
 723
 724             // something changed, copy seq, eval and return it
 725             epvector *s=new epvector;
 726             s->reserve(seq.size());
 727
 728             // copy parts of seq which are known not to have changed
 729             epvector::const_iterator cit2=seq.begin();
 730             while (cit2!=cit) {
 731                 s->push_back(*cit2);
 732                 ++cit2;
 733             }
 734             // copy first changed element
 735             s->push_back(split_ex_to_pair(expanded_factor));
 736             ++cit2;
 737             // copy rest
 738             while (cit2!=last) {
 739                 s->push_back(split_ex_to_pair(recombine_pair_to_ex(*cit2).expand(options)));
 740                 ++cit2;
 741             }
 742             return s;
 743         }
 744         ++cit;
 745     }
 746
 747     return 0; // nothing has changed
 748 }
 749
 750 //////////
 751 // static member variables
 752 //////////
 753
 754 // protected
 755
 756 unsigned mul::precedence=50;
 757
 758
 759 //////////
 760 // global constants
 761 //////////
 762
 763 const mul some_mul;
 764 type_info const & typeid_mul=typeid(some_mul);
 765
 766 #ifndef NO_GINAC_NAMESPACE
 767 } // namespace GiNaC
 768 #endif // ndef NO_GINAC_NAMESPACE