- changed old-style power() to new-style pow()
[ginac.git] / ginac / add.cpp
1 /** @file add.cpp
2  *
3  *  Implementation of GiNaC's sums 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 <iostream>
24 #include <stdexcept>
25
26 #include "add.h"
27 #include "mul.h"
28 #include "archive.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 GINAC_IMPLEMENT_REGISTERED_CLASS(add, expairseq)
37
38 //////////
39 // default constructor, destructor, copy constructor assignment operator and helpers
40 //////////
41
42 // public
43
44 add::add()
45 {
46     debugmsg("add default constructor",LOGLEVEL_CONSTRUCT);
47     tinfo_key = TINFO_add;
48 }
49
50 add::~add()
51 {
52     debugmsg("add destructor",LOGLEVEL_DESTRUCT);
53     destroy(0);
54 }
55
56 add::add(add const & other)
57 {
58     debugmsg("add copy constructor",LOGLEVEL_CONSTRUCT);
59     copy(other);
60 }
61
62 add const & add::operator=(add const & other)
63 {
64     debugmsg("add operator=",LOGLEVEL_ASSIGNMENT);
65     if (this != &other) {
66         destroy(1);
67         copy(other);
68     }
69     return *this;
70 }
71
72 // protected
73
74 void add::copy(add const & other)
75 {
76     inherited::copy(other);
77 }
78
79 void add::destroy(bool call_parent)
80 {
81     if (call_parent) inherited::destroy(call_parent);
82 }
83
84 //////////
85 // other constructors
86 //////////
87
88 // public
89
90 add::add(ex const & lh, ex const & rh)
91 {
92     debugmsg("add constructor from ex,ex",LOGLEVEL_CONSTRUCT);
93     tinfo_key = TINFO_add;
94     overall_coeff=_ex0();
95     construct_from_2_ex(lh,rh);
96     GINAC_ASSERT(is_canonical());
97 }
98
99 add::add(exvector const & v)
100 {
101     debugmsg("add constructor from exvector",LOGLEVEL_CONSTRUCT);
102     tinfo_key = TINFO_add;
103     overall_coeff=_ex0();
104     construct_from_exvector(v);
105     GINAC_ASSERT(is_canonical());
106 }
107
108 /*
109 add::add(epvector const & v, bool do_not_canonicalize)
110 {
111     debugmsg("add constructor from epvector,bool",LOGLEVEL_CONSTRUCT);
112     tinfo_key = TINFO_add;
113     if (do_not_canonicalize) {
114         seq=v;
115 #ifdef EXPAIRSEQ_USE_HASHTAB
116         combine_same_terms(); // to build hashtab
117 #endif // def EXPAIRSEQ_USE_HASHTAB
118     } else {
119         construct_from_epvector(v);
120     }
121     GINAC_ASSERT(is_canonical());
122 }
123 */
124
125 add::add(epvector const & v)
126 {
127     debugmsg("add constructor from epvector",LOGLEVEL_CONSTRUCT);
128     tinfo_key = TINFO_add;
129     overall_coeff=_ex0();
130     construct_from_epvector(v);
131     GINAC_ASSERT(is_canonical());
132 }
133
134 add::add(epvector const & v, ex const & oc)
135 {
136     debugmsg("add constructor from epvector,ex",LOGLEVEL_CONSTRUCT);
137     tinfo_key = TINFO_add;
138     overall_coeff=oc;
139     construct_from_epvector(v);
140     GINAC_ASSERT(is_canonical());
141 }
142
143 add::add(epvector * vp, ex const & oc)
144 {
145     debugmsg("add constructor from epvector *,ex",LOGLEVEL_CONSTRUCT);
146     tinfo_key = TINFO_add;
147     GINAC_ASSERT(vp!=0);
148     overall_coeff=oc;
149     construct_from_epvector(*vp);
150     delete vp;
151     GINAC_ASSERT(is_canonical());
152 }
153
154 //////////
155 // archiving
156 //////////
157
158 /** Construct object from archive_node. */
159 add::add(const archive_node &n, const lst &sym_lst) : inherited(n, sym_lst)
160 {
161     debugmsg("add constructor from archive_node", LOGLEVEL_CONSTRUCT);
162 }
163
164 /** Unarchive the object. */
165 ex add::unarchive(const archive_node &n, const lst &sym_lst)
166 {
167     return (new add(n, sym_lst))->setflag(status_flags::dynallocated);
168 }
169
170 /** Archive the object. */
171 void add::archive(archive_node &n) const
172 {
173     inherited::archive(n);
174 }
175
176 //////////
177 // functions overriding virtual functions from bases classes
178 //////////
179
180 // public
181
182 basic * add::duplicate() const
183 {
184     debugmsg("add duplicate",LOGLEVEL_DUPLICATE);
185     return new add(*this);
186 }
187
188 /*void add::print(ostream & os, unsigned upper_precedence) const
189 {
190     debugmsg("add print",LOGLEVEL_PRINT);
191     if (precedence<=upper_precedence) os << "(";
192     numeric coeff;
193     bool first=true;
194     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
195         coeff = ex_to_numeric(cit->coeff);
196         if (!first) {
197             if (coeff.csgn()==-1) os << '-'; else os << '+';
198         } else {
199             if (coeff.csgn()==-1) os << '-';
200             first=false;
201         }
202         if (!coeff.is_equal(_num1()) &&
203             !coeff.is_equal(_num_1())) {
204             if (coeff.is_rational()) {
205                 if (coeff.is_negative())
206                     os << -coeff;
207                 else
208                     os << coeff;
209             } else {
210                 if (coeff.csgn()==-1)
211                     (-coeff).print(os, precedence);
212                 else
213                     coeff.print(os, precedence);
214             }
215             os << '*';
216         }
217         os << cit->rest;
218     }
219     // print the overall numeric coefficient, if present:
220     if (!overall_coeff.is_zero()) {
221         if (overall_coeff.info(info_flags::positive)) os << '+';
222         os << overall_coeff;
223     }
224     if (precedence<=upper_precedence) os << ")";
225 }*/
226
227 void add::print(ostream & os, unsigned upper_precedence) const
228 {
229     debugmsg("add print",LOGLEVEL_PRINT);
230     if (precedence<=upper_precedence) os << "(";
231     numeric coeff;
232     bool first = true;
233     // first print the overall numeric coefficient, if present:
234     if (!overall_coeff.is_zero()) {
235         os << overall_coeff;
236         first = false;
237     }
238     // then proceed with the remaining factors:
239     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
240         coeff = ex_to_numeric(cit->coeff);
241         if (!first) {
242             if (coeff.csgn()==-1) os << '-'; else os << '+';
243         } else {
244             if (coeff.csgn()==-1) os << '-';
245             first = false;
246         }
247         if (!coeff.is_equal(_num1()) &&
248             !coeff.is_equal(_num_1())) {
249             if (coeff.is_rational()) {
250                 if (coeff.is_negative())
251                     os << -coeff;
252                 else
253                     os << coeff;
254             } else {
255                 if (coeff.csgn()==-1)
256                     (-coeff).print(os, precedence);
257                 else
258                     coeff.print(os, precedence);
259             }
260             os << '*';
261         }
262         os << cit->rest;
263     }
264     if (precedence<=upper_precedence) os << ")";
265 }
266
267 void add::printraw(ostream & os) const
268 {
269     debugmsg("add printraw",LOGLEVEL_PRINT);
270
271     os << "+(";
272     for (epvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) {
273         os << "(";
274         (*it).rest.bp->printraw(os);
275         os << ",";
276         (*it).coeff.bp->printraw(os);        
277         os << "),";
278     }
279     os << ",hash=" << hashvalue << ",flags=" << flags;
280     os << ")";
281 }
282
283 void add::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const
284 {
285     debugmsg("add print csrc", LOGLEVEL_PRINT);
286     if (precedence <= upper_precedence)
287         os << "(";
288
289     // Print arguments, separated by "+"
290     epvector::const_iterator it = seq.begin();
291     epvector::const_iterator itend = seq.end();
292     while (it != itend) {
293
294         // If the coefficient is -1, it is replaced by a single minus sign
295         if (it->coeff.compare(_num1()) == 0) {
296             it->rest.bp->printcsrc(os, type, precedence);
297         } else if (it->coeff.compare(_num_1()) == 0) {
298             os << "-";
299             it->rest.bp->printcsrc(os, type, precedence);
300         } else if (ex_to_numeric(it->coeff).numer().compare(_num1()) == 0) {
301             it->rest.bp->printcsrc(os, type, precedence);
302             os << "/";
303             ex_to_numeric(it->coeff).denom().printcsrc(os, type, precedence);
304         } else if (ex_to_numeric(it->coeff).numer().compare(_num_1()) == 0) {
305             os << "-";
306             it->rest.bp->printcsrc(os, type, precedence);
307             os << "/";
308             ex_to_numeric(it->coeff).denom().printcsrc(os, type, precedence);
309         } else {
310             it->coeff.bp->printcsrc(os, type, precedence);
311             os << "*";
312             it->rest.bp->printcsrc(os, type, precedence);
313         }
314
315         // Separator is "+", except if the following expression would have a leading minus sign
316         it++;
317         if (it != itend && !(it->coeff.compare(_num0()) < 0 || (it->coeff.compare(_num1()) == 0 && is_ex_exactly_of_type(it->rest, numeric) && it->rest.compare(_num0()) < 0)))
318             os << "+";
319     }
320     
321     if (!overall_coeff.is_equal(_ex0())) {
322         if (overall_coeff.info(info_flags::positive)) os << '+';
323         overall_coeff.bp->printcsrc(os,type,precedence);
324     }
325     
326     if (precedence <= upper_precedence)
327         os << ")";
328 }
329
330 bool add::info(unsigned inf) const
331 {
332     // TODO: optimize
333     if (inf==info_flags::polynomial ||
334         inf==info_flags::integer_polynomial ||
335         inf==info_flags::cinteger_polynomial ||
336         inf==info_flags::rational_polynomial ||
337         inf==info_flags::crational_polynomial ||
338         inf==info_flags::rational_function) {
339         for (epvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) {
340             if (!(recombine_pair_to_ex(*it).info(inf)))
341                 return false;
342         }
343         return overall_coeff.info(inf);
344     } else {
345         return inherited::info(inf);
346     }
347 }
348
349 int add::degree(symbol const & s) const
350 {
351     int deg=INT_MIN;
352     if (!overall_coeff.is_equal(_ex0())) {
353         deg=0;
354     }
355     int cur_deg;
356     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
357         cur_deg=(*cit).rest.degree(s);
358         if (cur_deg>deg) deg=cur_deg;
359     }
360     return deg;
361 }
362
363 int add::ldegree(symbol const & s) const
364 {
365     int deg=INT_MAX;
366     if (!overall_coeff.is_equal(_ex0())) {
367         deg=0;
368     }
369     int cur_deg;
370     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
371         cur_deg=(*cit).rest.ldegree(s);
372         if (cur_deg<deg) deg=cur_deg;
373     }
374     return deg;
375 }
376
377 ex add::coeff(symbol const & s, int const n) const
378 {
379     epvector coeffseq;
380     coeffseq.reserve(seq.size());
381
382     epvector::const_iterator it=seq.begin();
383     while (it!=seq.end()) {
384         coeffseq.push_back(combine_ex_with_coeff_to_pair((*it).rest.coeff(s,n),
385                                                          (*it).coeff));
386         ++it;
387     }
388     if (n==0) {
389         return (new add(coeffseq,overall_coeff))->setflag(status_flags::dynallocated);
390     }
391     return (new add(coeffseq))->setflag(status_flags::dynallocated);
392 }
393
394 ex add::eval(int level) const
395 {
396     // simplifications: +(;c) -> c
397     //                  +(x;1) -> x
398
399     debugmsg("add eval",LOGLEVEL_MEMBER_FUNCTION);
400
401     epvector * evaled_seqp=evalchildren(level);
402     if (evaled_seqp!=0) {
403         // do more evaluation later
404         return (new add(evaled_seqp,overall_coeff))->
405                    setflag(status_flags::dynallocated);
406     }
407
408 #ifdef DO_GINAC_ASSERT
409     for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
410         GINAC_ASSERT(!is_ex_exactly_of_type((*cit).rest,add));
411         if (is_ex_exactly_of_type((*cit).rest,numeric)) {
412             dbgprint();
413         }
414         GINAC_ASSERT(!is_ex_exactly_of_type((*cit).rest,numeric));
415     }
416 #endif // def DO_GINAC_ASSERT
417
418     if (flags & status_flags::evaluated) {
419         GINAC_ASSERT(seq.size()>0);
420         GINAC_ASSERT((seq.size()>1)||!overall_coeff.is_equal(_ex0()));
421         return *this;
422     }
423
424     int seq_size=seq.size();
425     if (seq_size==0) {
426         // +(;c) -> c
427         return overall_coeff;
428     } else if ((seq_size==1)&&overall_coeff.is_equal(_ex0())) {
429         // +(x;0) -> x
430         return recombine_pair_to_ex(*(seq.begin()));
431     }
432     return this->hold();
433 }
434
435 exvector add::get_indices(void) const
436 {
437     // FIXME: all terms in the sum should have the same indices (compatible
438     // tensors) however this is not checked, since there is no function yet
439     // which compares indices (idxvector can be unsorted)
440     if (seq.size()==0) {
441         return exvector();
442     }
443     return (seq.begin())->rest.get_indices();
444 }    
445
446 ex add::simplify_ncmul(exvector const & v) const
447 {
448     if (seq.size()==0) {
449         return inherited::simplify_ncmul(v);
450     }
451     return (*seq.begin()).rest.simplify_ncmul(v);
452 }    
453
454 // protected
455
456 int add::compare_same_type(basic const & other) const
457 {
458     return inherited::compare_same_type(other);
459 }
460
461 bool add::is_equal_same_type(basic const & other) const
462 {
463     return inherited::is_equal_same_type(other);
464 }
465
466 unsigned add::return_type(void) const
467 {
468     if (seq.size()==0) {
469         return return_types::commutative;
470     }
471     return (*seq.begin()).rest.return_type();
472 }
473    
474 unsigned add::return_type_tinfo(void) const
475 {
476     if (seq.size()==0) {
477         return tinfo_key;
478     }
479     return (*seq.begin()).rest.return_type_tinfo();
480 }
481
482 ex add::thisexpairseq(epvector const & v, ex const & oc) const
483 {
484     return (new add(v,oc))->setflag(status_flags::dynallocated);
485 }
486
487 ex add::thisexpairseq(epvector * vp, ex const & oc) const
488 {
489     return (new add(vp,oc))->setflag(status_flags::dynallocated);
490 }
491
492 expair add::split_ex_to_pair(ex const & e) const
493 {
494     if (is_ex_exactly_of_type(e,mul)) {
495         mul const & mulref=ex_to_mul(e);
496         ex numfactor=mulref.overall_coeff;
497         // mul * mulcopyp=static_cast<mul *>(mulref.duplicate());
498         mul * mulcopyp=new mul(mulref);
499         mulcopyp->overall_coeff=_ex1();
500         mulcopyp->clearflag(status_flags::evaluated);
501         mulcopyp->clearflag(status_flags::hash_calculated);
502         return expair(mulcopyp->setflag(status_flags::dynallocated),numfactor);
503     }
504     return expair(e,_ex1());
505 }
506
507 expair add::combine_ex_with_coeff_to_pair(ex const & e,
508                                           ex const & c) const
509 {
510     GINAC_ASSERT(is_ex_exactly_of_type(c,numeric));
511     if (is_ex_exactly_of_type(e,mul)) {
512         mul const & mulref=ex_to_mul(e);
513         ex numfactor=mulref.overall_coeff;
514         //mul * mulcopyp=static_cast<mul *>(mulref.duplicate());
515         mul * mulcopyp=new mul(mulref);
516         mulcopyp->overall_coeff=_ex1();
517         mulcopyp->clearflag(status_flags::evaluated);
518         mulcopyp->clearflag(status_flags::hash_calculated);
519         if (are_ex_trivially_equal(c,_ex1())) {
520             return expair(mulcopyp->setflag(status_flags::dynallocated),numfactor);
521         } else if (are_ex_trivially_equal(numfactor,_ex1())) {
522             return expair(mulcopyp->setflag(status_flags::dynallocated),c);
523         }
524         return expair(mulcopyp->setflag(status_flags::dynallocated),
525                           ex_to_numeric(numfactor).mul_dyn(ex_to_numeric(c)));
526     } else if (is_ex_exactly_of_type(e,numeric)) {
527         if (are_ex_trivially_equal(c,_ex1())) {
528             return expair(e,_ex1());
529         }
530         return expair(ex_to_numeric(e).mul_dyn(ex_to_numeric(c)),_ex1());
531     }
532     return expair(e,c);
533 }
534     
535 expair add::combine_pair_with_coeff_to_pair(expair const & p,
536                                             ex const & c) const
537 {
538     GINAC_ASSERT(is_ex_exactly_of_type(p.coeff,numeric));
539     GINAC_ASSERT(is_ex_exactly_of_type(c,numeric));
540
541     if (is_ex_exactly_of_type(p.rest,numeric)) {
542         GINAC_ASSERT(ex_to_numeric(p.coeff).is_equal(_num1())); // should be normalized
543         return expair(ex_to_numeric(p.rest).mul_dyn(ex_to_numeric(c)),_ex1());
544     }
545
546     return expair(p.rest,ex_to_numeric(p.coeff).mul_dyn(ex_to_numeric(c)));
547 }
548     
549 ex add::recombine_pair_to_ex(expair const & p) const
550 {
551     //if (p.coeff.compare(_ex1())==0) {
552     //if (are_ex_trivially_equal(p.coeff,_ex1())) {
553     if (ex_to_numeric(p.coeff).is_equal(_num1())) {
554         return p.rest;
555     } else {
556         return p.rest*p.coeff;
557     }
558 }
559
560 ex add::expand(unsigned options) const
561 {
562     epvector * vp=expandchildren(options);
563     if (vp==0) {
564         return *this;
565     }
566     return (new add(vp,overall_coeff))->setflag(status_flags::expanded    |
567                                                 status_flags::dynallocated );
568 }
569
570 //////////
571 // new virtual functions which can be overridden by derived classes
572 //////////
573
574 // none
575
576 //////////
577 // non-virtual functions in this class
578 //////////
579
580 // none
581
582 //////////
583 // static member variables
584 //////////
585
586 // protected
587
588 unsigned add::precedence=40;
589
590 //////////
591 // global constants
592 //////////
593
594 const add some_add;
595 type_info const & typeid_add=typeid(some_add);
596
597 #ifndef NO_GINAC_NAMESPACE
598 } // namespace GiNaC
599 #endif // ndef NO_GINAC_NAMESPACE