1 /** @file basic.cpp
2 *
3 * Implementation of GiNaC's ABC. */
4
5 /*
6 * GiNaC Copyright (C) 1999-2022 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., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
21 */
22
23 #include "basic.h"
24 #include "ex.h"
25 #include "numeric.h"
26 #include "power.h"
27 #include "add.h"
28 #include "symbol.h"
29 #include "lst.h"
30 #include "ncmul.h"
31 #include "relational.h"
32 #include "operators.h"
33 #include "wildcard.h"
34 #include "archive.h"
35 #include "utils.h"
36 #include "hash_seed.h"
37 #include "inifcns.h"
38
39 #include <iostream>
40 #include <stdexcept>
41 #include <typeinfo>
42
43 namespace GiNaC {
44
45 GINAC_IMPLEMENT_REGISTERED_CLASS_OPT(basic, void,
46 print_func<print_context>(&basic::do_print).
47 print_func<print_tree>(&basic::do_print_tree).
48 print_func<print_python_repr>(&basic::do_print_python_repr))
49
50 //////////
51 // default constructor, destructor, copy constructor and assignment operator
52 //////////
53
54 // public
55
56 /** basic copy constructor: implicitly assumes that the other class is of
57 * the exact same type (as it's used by duplicate()), so it can copy the
58 * tinfo_key and the hash value. */
basic(const basic & other)59 basic::basic(const basic & other) : flags(other.flags & ~status_flags::dynallocated), hashvalue(other.hashvalue)
60 {
61 }
62
63 /** basic assignment operator: the other object might be of a derived class. */
operator =(const basic & other)64 const basic & basic::operator=(const basic & other)
65 {
66 unsigned fl = other.flags & ~status_flags::dynallocated;
67 if (typeid(*this) != typeid(other)) {
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 fl &= ~(status_flags::evaluated | status_flags::expanded | status_flags::hash_calculated);
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 set_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. */
read_archive(const archive_node & n,lst & syms)96 void basic::read_archive(const archive_node& n, lst& syms)
97 { }
98
99 /** Archive the object. */
archive(archive_node & n) const100 void basic::archive(archive_node &n) const
101 {
102 n.add_string("class", class_name());
103 }
104
105 //////////
106 // new virtual functions which can be overridden by derived classes
107 //////////
108
109 // public
110
111 /** Output to stream. This performs double dispatch on the dynamic type of
112 * *this and the dynamic type of the supplied print context.
113 * @param c print context object that describes the output formatting
114 * @param level value that is used to identify the precedence or indentation
115 * level for placing parentheses and formatting */
print(const print_context & c,unsigned level) const116 void basic::print(const print_context & c, unsigned level) const
117 {
118 print_dispatch(get_class_info(), c, level);
119 }
120
121 /** Like print(), but dispatch to the specified class. Can be used by
122 * implementations of print methods to dispatch to the method of the
123 * superclass.
124 *
125 * @see basic::print */
print_dispatch(const registered_class_info & ri,const print_context & c,unsigned level) const126 void basic::print_dispatch(const registered_class_info & ri, 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 = &ri;
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. */
do_print(const print_context & c,unsigned level) const169 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. */
do_print_tree(const print_tree & c,unsigned level) const175 void basic::do_print_tree(const print_tree & c, unsigned level) const
176 {
177 c.s << std::string(level, ' ') << class_name() << " @" << this
178 << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec;
179 if (nops())
180 c.s << ", nops=" << nops();
181 c.s << std::endl;
182 for (size_t i=0; i<nops(); ++i)
183 op(i).print(c, level + c.delta_indent);
184 }
185
186 /** Python parsable output to stream. */
do_print_python_repr(const print_python_repr & c,unsigned level) const187 void basic::do_print_python_repr(const print_python_repr & c, unsigned level) const
188 {
189 c.s << class_name() << "()";
190 }
191
192 /** Little wrapper around print to be called within a debugger.
193 * This is needed because you cannot call foo.print(cout) from within the
194 * debugger because it might not know what cout is. This method can be
195 * invoked with no argument and it will simply print to stdout.
196 *
197 * @see basic::print
198 * @see basic::dbgprinttree */
dbgprint() const199 void basic::dbgprint() const
200 {
201 this->print(print_dflt(std::cerr));
202 std::cerr << std::endl;
203 }
204
205 /** Little wrapper around printtree to be called within a debugger.
206 *
207 * @see basic::dbgprint */
dbgprinttree() const208 void basic::dbgprinttree() const
209 {
210 this->print(print_tree(std::cerr));
211 }
212
213 /** Return relative operator precedence (for parenthezing output). */
precedence() const214 unsigned basic::precedence() const
215 {
216 return 70;
217 }
218
219 /** Information about the object.
220 *
221 * @see class info_flags */
info(unsigned inf) const222 bool basic::info(unsigned inf) const
223 {
224 // all possible properties are false for basic objects
225 return false;
226 }
227
228 /** Number of operands/members. */
nops() const229 size_t basic::nops() const
230 {
231 // iterating from 0 to nops() on atomic objects should be an empty loop,
232 // and accessing their elements is a range error. Container objects should
233 // override this.
234 return 0;
235 }
236
237 /** Return operand/member at position i. */
op(size_t i) const238 ex basic::op(size_t i) const
239 {
240 throw(std::range_error(std::string("basic::op(): ") + class_name() + std::string(" has no operands")));
241 }
242
243 /** Return modifiable operand/member at position i. */
let_op(size_t i)244 ex & basic::let_op(size_t i)
245 {
246 ensure_if_modifiable();
247 throw(std::range_error(std::string("basic::let_op(): ") + class_name() + std::string(" has no operands")));
248 }
249
operator [](const ex & index) const250 ex basic::operator[](const ex & index) const
251 {
252 if (is_exactly_a<numeric>(index))
253 return op(static_cast<size_t>(ex_to<numeric>(index).to_int()));
254
255 throw(std::invalid_argument(std::string("non-numeric indices not supported by ") + class_name()));
256 }
257
operator [](size_t i) const258 ex basic::operator[](size_t i) const
259 {
260 return op(i);
261 }
262
operator [](const ex & index)263 ex & basic::operator[](const ex & index)
264 {
265 if (is_exactly_a<numeric>(index))
266 return let_op(ex_to<numeric>(index).to_int());
267
268 throw(std::invalid_argument(std::string("non-numeric indices not supported by ") + class_name()));
269 }
270
operator [](size_t i)271 ex & basic::operator[](size_t i)
272 {
273 return let_op(i);
274 }
275
276 /** Test for occurrence of a pattern. An object 'has' a pattern if it matches
277 * the pattern itself or one of the children 'has' it. As a consequence
278 * (according to the definition of children) given e=x+y+z, e.has(x) is true
279 * but e.has(x+y) is false. */
has(const ex & pattern,unsigned options) const280 bool basic::has(const ex & pattern, unsigned options) const
281 {
282 exmap repl_lst;
283 if (match(pattern, repl_lst))
284 return true;
285 for (size_t i=0; i<nops(); i++)
286 if (op(i).has(pattern, options))
287 return true;
288
289 return false;
290 }
291
292 /** Construct new expression by applying the specified function to all
293 * sub-expressions (one level only, not recursively). */
map(map_function & f) const294 ex basic::map(map_function & f) const
295 {
296 size_t num = nops();
297 if (num == 0)
298 return *this;
299
300 basic *copy = nullptr;
301 for (size_t i=0; i<num; i++) {
302 const ex & o = op(i);
303 const ex & n = f(o);
304 if (!are_ex_trivially_equal(o, n)) {
305 if (copy == nullptr)
306 copy = duplicate();
307 copy->let_op(i) = n;
308 }
309 }
310
311 if (copy) {
312 copy->clearflag(status_flags::hash_calculated | status_flags::expanded);
313 return *copy;
314 } else
315 return *this;
316 }
317
318 /** Check whether this is a polynomial in the given variables. */
is_polynomial(const ex & var) const319 bool basic::is_polynomial(const ex & var) const
320 {
321 return !has(var) || is_equal(ex_to<basic>(var));
322 }
323
324 /** Return degree of highest power in object s. */
degree(const ex & s) const325 int basic::degree(const ex & s) const
326 {
327 return is_equal(ex_to<basic>(s)) ? 1 : 0;
328 }
329
330 /** Return degree of lowest power in object s. */
ldegree(const ex & s) const331 int basic::ldegree(const ex & s) const
332 {
333 return is_equal(ex_to<basic>(s)) ? 1 : 0;
334 }
335
336 /** Return coefficient of degree n in object s. */
coeff(const ex & s,int n) const337 ex basic::coeff(const ex & s, int n) const
338 {
339 if (is_equal(ex_to<basic>(s)))
340 return n==1 ? _ex1 : _ex0;
341 else
342 return n==0 ? *this : _ex0;
343 }
344
345 /** Sort expanded expression in terms of powers of some object(s).
346 * @param s object(s) to sort in
347 * @param distributed recursive or distributed form (only used when s is a list) */
collect(const ex & s,bool distributed) const348 ex basic::collect(const ex & s, bool distributed) const
349 {
350 ex x;
351 if (is_a<lst>(s)) {
352
353 // List of objects specified
354 if (s.nops() == 0)
355 return *this;
356 if (s.nops() == 1)
357 return collect(s.op(0));
358
359 else if (distributed) {
360
361 x = this->expand();
362 if (! is_a<add>(x))
363 return x;
364 const lst& l(ex_to<lst>(s));
365
366 exmap cmap;
367 cmap[_ex1] = _ex0;
368 for (const auto & xi : x) {
369 ex key = _ex1;
370 ex pre_coeff = xi;
371 for (auto & li : l) {
372 int cexp = pre_coeff.degree(li);
373 pre_coeff = pre_coeff.coeff(li, cexp);
374 key *= pow(li, cexp);
375 }
376 auto ci = cmap.find(key);
377 if (ci != cmap.end())
378 ci->second += pre_coeff;
379 else
380 cmap.insert(exmap::value_type(key, pre_coeff));
381 }
382
383 exvector resv;
384 for (auto & mi : cmap)
385 resv.push_back((mi.first)*(mi.second));
386 return dynallocate<add>(resv);
387
388 } else {
389
390 // Recursive form
391 x = *this;
392 size_t n = s.nops() - 1;
393 while (true) {
394 x = x.collect(s[n]);
395 if (n == 0)
396 break;
397 n--;
398 }
399 }
400
401 } else {
402
403 // Only one object specified
404 for (int n=this->ldegree(s); n<=this->degree(s); ++n)
405 x += this->coeff(s,n)*power(s,n);
406 }
407
408 // correct for lost fractional arguments and return
409 return x + (*this - x).expand();
410 }
411
412 /** Perform automatic non-interruptive term rewriting rules. */
eval() const413 ex basic::eval() const
414 {
415 // There is nothing to do for basic objects:
416 return hold();
417 }
418
419 /** Function object to be applied by basic::evalf(). */
420 struct evalf_map_function : public map_function {
operator ()GiNaC::evalf_map_function421 ex operator()(const ex & e) override { return evalf(e); }
422 };
423
424 /** Evaluate object numerically. */
evalf() const425 ex basic::evalf() const
426 {
427 if (nops() == 0)
428 return *this;
429 else {
430 evalf_map_function map_evalf;
431 return map(map_evalf);
432 }
433 }
434
435 /** Function object to be applied by basic::evalm(). */
436 struct evalm_map_function : public map_function {
operator ()GiNaC::evalm_map_function437 ex operator()(const ex & e) override { return evalm(e); }
438 } map_evalm;
439
440 /** Evaluate sums, products and integer powers of matrices. */
evalm() const441 ex basic::evalm() const
442 {
443 if (nops() == 0)
444 return *this;
445 else
446 return map(map_evalm);
447 }
448
449 /** Function object to be applied by basic::eval_integ(). */
450 struct eval_integ_map_function : public map_function {
operator ()GiNaC::eval_integ_map_function451 ex operator()(const ex & e) override { return eval_integ(e); }
452 } map_eval_integ;
453
454 /** Evaluate integrals, if result is known. */
eval_integ() const455 ex basic::eval_integ() const
456 {
457 if (nops() == 0)
458 return *this;
459 else
460 return map(map_eval_integ);
461 }
462
463 /** Perform automatic symbolic evaluations on indexed expression that
464 * contains this object as the base expression. */
eval_indexed(const basic & i) const465 ex basic::eval_indexed(const basic & i) const
466 // this function can't take a "const ex & i" because that would result
467 // in an infinite eval() loop
468 {
469 // There is nothing to do for basic objects
470 return i.hold();
471 }
472
473 /** Add two indexed expressions. They are guaranteed to be of class indexed
474 * (or a subclass) and their indices are compatible. This function is used
475 * internally by simplify_indexed().
476 *
477 * @param self First indexed expression; its base object is *this
478 * @param other Second indexed expression
479 * @return sum of self and other
480 * @see ex::simplify_indexed() */
add_indexed(const ex & self,const ex & other) const481 ex basic::add_indexed(const ex & self, const ex & other) const
482 {
483 return self + other;
484 }
485
486 /** Multiply an indexed expression with a scalar. This function is used
487 * internally by simplify_indexed().
488 *
489 * @param self Indexed expression; its base object is *this
490 * @param other Numeric value
491 * @return product of self and other
492 * @see ex::simplify_indexed() */
scalar_mul_indexed(const ex & self,const numeric & other) const493 ex basic::scalar_mul_indexed(const ex & self, const numeric & other) const
494 {
495 return self * other;
496 }
497
498 /** Try to contract two indexed expressions that appear in the same product.
499 * If a contraction exists, the function overwrites one or both of the
500 * expressions and returns true. Otherwise it returns false. It is
501 * guaranteed that both expressions are of class indexed (or a subclass)
502 * and that at least one dummy index has been found. This functions is
503 * used internally by simplify_indexed().
504 *
505 * @param self Pointer to first indexed expression; its base object is *this
506 * @param other Pointer to second indexed expression
507 * @param v The complete vector of factors
508 * @return true if the contraction was successful, false otherwise
509 * @see ex::simplify_indexed() */
contract_with(exvector::iterator self,exvector::iterator other,exvector & v) const510 bool basic::contract_with(exvector::iterator self, exvector::iterator other, exvector & v) const
511 {
512 // Do nothing
513 return false;
514 }
515
516 /** Check whether the expression matches a given pattern. For every wildcard
517 * object in the pattern, a pair with the wildcard as a key and matching
518 * expression as a value is added to repl_lst. */
match(const ex & pattern,exmap & repl_lst) const519 bool basic::match(const ex & pattern, exmap& repl_lst) const
520 {
521 /*
522 Sweet sweet shapes, sweet sweet shapes,
523 That's the key thing, right right.
524 Feed feed face, feed feed shapes,
525 But who is the king tonight?
526 Who is the king tonight?
527 Pattern is the thing, the key thing-a-ling,
528 But who is the king of Pattern?
529 But who is the king, the king thing-a-ling,
530 Who is the king of Pattern?
531 Bog is the king, the king thing-a-ling,
532 Bog is the king of Pattern.
533 Ba bu-bu-bu-bu bu-bu-bu-bu-bu-bu bu-bu
534 Bog is the king of Pattern.
535 */
536
537 if (is_exactly_a<wildcard>(pattern)) {
538
539 // Wildcard matches anything, but check whether we already have found
540 // a match for that wildcard first (if so, the earlier match must be
541 // the same expression)
542 for (auto & it : repl_lst) {
543 if (it.first.is_equal(pattern))
544 return is_equal(ex_to<basic>(it.second));
545 }
546 repl_lst[pattern] = *this;
547 return true;
548
549 } else {
550
551 // Expression must be of the same type as the pattern
552 if (typeid(*this) != typeid(ex_to<basic>(pattern)))
553 return false;
554
555 // Number of subexpressions must match
556 if (nops() != pattern.nops())
557 return false;
558
559 // No subexpressions? Then just compare the objects (there can't be
560 // wildcards in the pattern)
561 if (nops() == 0)
562 return is_equal_same_type(ex_to<basic>(pattern));
563
564 // Check whether attributes that are not subexpressions match
565 if (!match_same_type(ex_to<basic>(pattern)))
566 return false;
567
568 // Even if the expression does not match the pattern, some of
569 // its subexpressions could match it. For example, x^5*y^(-1)
570 // does not match the pattern $0^5, but its subexpression x^5
571 // does. So, save repl_lst in order to not add bogus entries.
572 exmap tmp_repl = repl_lst;
573 // Otherwise the subexpressions must match one-to-one
574 for (size_t i=0; i<nops(); i++)
575 if (!op(i).match(pattern.op(i), tmp_repl))
576 return false;
577
578 // Looks similar enough, match found
579 repl_lst = tmp_repl;
580 return true;
581 }
582 }
583
584 /** Helper function for subs(). Does not recurse into subexpressions. */
subs_one_level(const exmap & m,unsigned options) const585 ex basic::subs_one_level(const exmap & m, unsigned options) const
586 {
587 if (options & subs_options::no_pattern) {
588 ex thisex = *this; // NB: *this may be deleted here.
589 auto it = m.find(thisex);
590 if (it != m.end())
591 return it->second;
592 return thisex;
593 } else {
594 for (auto & it : m) {
595 exmap repl_lst;
596 if (match(ex_to<basic>(it.first), repl_lst))
597 return it.second.subs(repl_lst, options | subs_options::no_pattern);
598 // avoid infinite recursion when re-substituting the wildcards
599 }
600 }
601
602 return *this;
603 }
604
605 /** Substitute a set of objects by arbitrary expressions. The ex returned
606 * will already be evaluated. */
subs(const exmap & m,unsigned options) const607 ex basic::subs(const exmap & m, unsigned options) const
608 {
609 size_t num = nops();
610 if (num) {
611
612 // Substitute in subexpressions
613 for (size_t i=0; i<num; i++) {
614 const ex & orig_op = op(i);
615 const ex & subsed_op = orig_op.subs(m, options);
616 if (!are_ex_trivially_equal(orig_op, subsed_op)) {
617
618 // Something changed, clone the object
619 basic *copy = duplicate();
620 copy->clearflag(status_flags::hash_calculated | status_flags::expanded);
621
622 // Substitute the changed operand
623 copy->let_op(i++) = subsed_op;
624
625 // Substitute the other operands
626 for (; i<num; i++)
627 copy->let_op(i) = op(i).subs(m, options);
628
629 // Perform substitutions on the new object as a whole
630 return copy->subs_one_level(m, options);
631 }
632 }
633 }
634
635 // Nothing changed or no subexpressions
636 return subs_one_level(m, options);
637 }
638
639 /** Default interface of nth derivative ex::diff(s, n). It should be called
640 * instead of ::derivative(s) for first derivatives and for nth derivatives it
641 * just recurses down.
642 *
643 * @param s symbol to differentiate in
644 * @param nth order of differentiation
645 * @see ex::diff */
diff(const symbol & s,unsigned nth) const646 ex basic::diff(const symbol & s, unsigned nth) const
647 {
648 // trivial: zeroth derivative
649 if (nth==0)
650 return ex(*this);
651
652 // evaluate unevaluated *this before differentiating
653 if (!(flags & status_flags::evaluated))
654 return ex(*this).diff(s, nth);
655
656 ex ndiff = this->derivative(s);
657 while (!ndiff.is_zero() && // stop differentiating zeros
658 nth>1) {
659 ndiff = ndiff.diff(s);
660 --nth;
661 }
662 return ndiff;
663 }
664
665 /** Return a vector containing the free indices of an expression. */
get_free_indices() const666 exvector basic::get_free_indices() const
667 {
668 return exvector(); // return an empty exvector
669 }
670
conjugate() const671 ex basic::conjugate() const
672 {
673 return *this;
674 }
675
real_part() const676 ex basic::real_part() const
677 {
678 return real_part_function(*this).hold();
679 }
680
imag_part() const681 ex basic::imag_part() const
682 {
683 return imag_part_function(*this).hold();
684 }
685
eval_ncmul(const exvector & v) const686 ex basic::eval_ncmul(const exvector & v) const
687 {
688 return hold_ncmul(v);
689 }
690
691 // protected
692
693 /** Function object to be applied by basic::derivative(). */
694 struct derivative_map_function : public map_function {
695 const symbol &s;
derivative_map_functionGiNaC::derivative_map_function696 derivative_map_function(const symbol &sym) : s(sym) {}
operator ()GiNaC::derivative_map_function697 ex operator()(const ex & e) override { return diff(e, s); }
698 };
699
700 /** Default implementation of ex::diff(). It maps the operation on the
701 * operands (or returns 0 when the object has no operands).
702 *
703 * @see ex::diff */
derivative(const symbol & s) const704 ex basic::derivative(const symbol & s) const
705 {
706 if (nops() == 0)
707 return _ex0;
708 else {
709 derivative_map_function map_derivative(s);
710 return map(map_derivative);
711 }
712 }
713
714 /** Returns order relation between two objects of same type. This needs to be
715 * implemented by each class. It may never return anything else than 0,
716 * signalling equality, or +1 and -1 signalling inequality and determining
717 * the canonical ordering. (Perl hackers will wonder why C++ doesn't feature
718 * the spaceship operator <=> for denoting just this.) */
compare_same_type(const basic & other) const719 int basic::compare_same_type(const basic & other) const
720 {
721 return compare_pointers(this, &other);
722 }
723
724 /** Returns true if two objects of same type are equal. Normally needs
725 * not be reimplemented as long as it wasn't overwritten by some parent
726 * class, since it just calls compare_same_type(). The reason why this
727 * function exists is that sometimes it is easier to determine equality
728 * than an order relation and then it can be overridden. */
is_equal_same_type(const basic & other) const729 bool basic::is_equal_same_type(const basic & other) const
730 {
731 return compare_same_type(other)==0;
732 }
733
734 /** Returns true if the attributes of two objects are similar enough for
735 * a match. This function must not match subexpressions (this is already
736 * done by basic::match()). Only attributes not accessible by op() should
737 * be compared. This is also the reason why this function doesn't take the
738 * wildcard replacement list from match() as an argument: only subexpressions
739 * are subject to wildcard matches. Also, this function only needs to be
740 * implemented for container classes because is_equal_same_type() is
741 * automatically used instead of match_same_type() if nops() == 0.
742 *
743 * @see basic::match */
match_same_type(const basic & other) const744 bool basic::match_same_type(const basic & other) const
745 {
746 // The default is to only consider subexpressions, but not any other
747 // attributes
748 return true;
749 }
750
return_type() const751 unsigned basic::return_type() const
752 {
753 return return_types::commutative;
754 }
755
return_type_tinfo() const756 return_type_t basic::return_type_tinfo() const
757 {
758 return_type_t rt;
759 rt.tinfo = &typeid(*this);
760 rt.rl = 0;
761 return rt;
762 }
763
764 /** Compute the hash value of an object and if it makes sense to store it in
765 * the objects status_flags, do so. The method inherited from class basic
766 * computes a hash value based on the type and hash values of possible
767 * members. For this reason it is well suited for container classes but
768 * atomic classes should override this implementation because otherwise they
769 * would all end up with the same hashvalue. */
calchash() const770 unsigned basic::calchash() const
771 {
772 unsigned v = make_hash_seed(typeid(*this));
773 for (size_t i=0; i<nops(); i++) {
774 v = rotate_left(v);
775 v ^= this->op(i).gethash();
776 }
777
778 // store calculated hash value only if object is already evaluated
779 if (flags & status_flags::evaluated) {
780 setflag(status_flags::hash_calculated);
781 hashvalue = v;
782 }
783
784 return v;
785 }
786
787 /** Function object to be applied by basic::expand(). */
788 struct expand_map_function : public map_function {
789 unsigned options;
expand_map_functionGiNaC::expand_map_function790 expand_map_function(unsigned o) : options(o) {}
operator ()GiNaC::expand_map_function791 ex operator()(const ex & e) override { return e.expand(options); }
792 };
793
794 /** Expand expression, i.e. multiply it out and return the result as a new
795 * expression. */
expand(unsigned options) const796 ex basic::expand(unsigned options) const
797 {
798 if (nops() == 0)
799 return (options == 0) ? setflag(status_flags::expanded) : *this;
800 else {
801 expand_map_function map_expand(options);
802 return ex_to<basic>(map(map_expand)).setflag(options == 0 ? status_flags::expanded : 0);
803 }
804 }
805
806
807 //////////
808 // non-virtual functions in this class
809 //////////
810
811 // public
812
813 /** Compare objects syntactically to establish canonical ordering.
814 * All compare functions return: -1 for *this less than other, 0 equal,
815 * 1 greater. */
compare(const basic & other) const816 int basic::compare(const basic & other) const
817 {
818 #ifdef GINAC_COMPARE_STATISTICS
819 compare_statistics.total_basic_compares++;
820 #endif
821 const unsigned hash_this = gethash();
822 const unsigned hash_other = other.gethash();
823 if (hash_this<hash_other) return -1;
824 if (hash_this>hash_other) return 1;
825 #ifdef GINAC_COMPARE_STATISTICS
826 compare_statistics.compare_same_hashvalue++;
827 #endif
828
829 const std::type_info& typeid_this = typeid(*this);
830 const std::type_info& typeid_other = typeid(other);
831 if (typeid_this == typeid_other) {
832 // int cmpval = compare_same_type(other);
833 // if (cmpval!=0) {
834 // std::cout << "hash collision, same type: "
835 // << *this << " and " << other << std::endl;
836 // this->print(print_tree(std::cout));
837 // std::cout << " and ";
838 // other.print(print_tree(std::cout));
839 // std::cout << std::endl;
840 // }
841 // return cmpval;
842 #ifdef GINAC_COMPARE_STATISTICS
843 compare_statistics.compare_same_type++;
844 #endif
845 return compare_same_type(other);
846 } else {
847 // std::cout << "hash collision, different types: "
848 // << *this << " and " << other << std::endl;
849 // this->print(print_tree(std::cout));
850 // std::cout << " and ";
851 // other.print(print_tree(std::cout));
852 // std::cout << std::endl;
853 return (typeid_this.before(typeid_other) ? -1 : 1);
854 }
855 }
856
857 /** Test for syntactic equality.
858 * This is only a quick test, meaning objects should be in the same domain.
859 * You might have to .expand(), .normal() objects first, depending on the
860 * domain of your computation, to get a more reliable answer.
861 *
862 * @see is_equal_same_type */
is_equal(const basic & other) const863 bool basic::is_equal(const basic & other) const
864 {
865 #ifdef GINAC_COMPARE_STATISTICS
866 compare_statistics.total_basic_is_equals++;
867 #endif
868 if (this->gethash()!=other.gethash())
869 return false;
870 #ifdef GINAC_COMPARE_STATISTICS
871 compare_statistics.is_equal_same_hashvalue++;
872 #endif
873 if (typeid(*this) != typeid(other))
874 return false;
875
876 #ifdef GINAC_COMPARE_STATISTICS
877 compare_statistics.is_equal_same_type++;
878 #endif
879 return is_equal_same_type(other);
880 }
881
882 // protected
883
884 /** Stop further evaluation.
885 *
886 * @see basic::eval */
hold() const887 const basic & basic::hold() const
888 {
889 return setflag(status_flags::evaluated);
890 }
891
892 /** Ensure the object may be modified without hurting others, throws if this
893 * is not the case. */
ensure_if_modifiable() const894 void basic::ensure_if_modifiable() const
895 {
896 if (get_refcount() > 1)
897 throw(std::runtime_error("cannot modify multiply referenced object"));
898 clearflag(status_flags::hash_calculated | status_flags::evaluated);
899 }
900
901 //////////
902 // global variables
903 //////////
904
905 #ifdef GINAC_COMPARE_STATISTICS
~compare_statistics_t()906 compare_statistics_t::~compare_statistics_t()
907 {
908 std::clog << "ex::compare() called " << total_compares << " times" << std::endl;
909 std::clog << "nontrivial compares: " << nontrivial_compares << " times" << std::endl;
910 std::clog << "basic::compare() called " << total_basic_compares << " times" << std::endl;
911 std::clog << "same hashvalue in compare(): " << compare_same_hashvalue << " times" << std::endl;
912 std::clog << "compare_same_type() called " << compare_same_type << " times" << std::endl;
913 std::clog << std::endl;
914 std::clog << "ex::is_equal() called " << total_is_equals << " times" << std::endl;
915 std::clog << "nontrivial is_equals: " << nontrivial_is_equals << " times" << std::endl;
916 std::clog << "basic::is_equal() called " << total_basic_is_equals << " times" << std::endl;
917 std::clog << "same hashvalue in is_equal(): " << is_equal_same_hashvalue << " times" << std::endl;
918 std::clog << "is_equal_same_type() called " << is_equal_same_type << " times" << std::endl;
919 std::clog << std::endl;
920 std::clog << "basic::gethash() called " << total_gethash << " times" << std::endl;
921 std::clog << "used cached hashvalue " << gethash_cached << " times" << std::endl;
922 }
923
924 compare_statistics_t compare_statistics;
925 #endif
926
927 } // namespace GiNaC
928