GiNaC  1.8.0
expairseq.cpp
Go to the documentation of this file.
1 
5 /*
6  * GiNaC Copyright (C) 1999-2020 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 "expairseq.h"
24 #include "lst.h"
25 #include "add.h"
26 #include "mul.h"
27 #include "power.h"
28 #include "relational.h"
29 #include "wildcard.h"
30 #include "archive.h"
31 #include "operators.h"
32 #include "utils.h"
33 #include "hash_seed.h"
34 #include "indexed.h"
35 #include "compiler.h"
36 
37 #include <algorithm>
38 #include <iostream>
39 #include <iterator>
40 #include <memory>
41 #include <stdexcept>
42 #include <string>
43 
44 namespace GiNaC {
45 
46 
49  print_func<print_tree>(&expairseq::do_print_tree))
50 
51 
52 // helper classes
55 
56 class epp_is_less
57 {
58 public:
59  bool operator()(const epp &lh, const epp &rh) const
60  {
61  return (*lh).is_less(*rh);
62  }
63 };
64 
66 // default constructor
68 
69 // public
70 
72 {}
73 
74 // protected
75 
77 // other constructors
79 
80 expairseq::expairseq(const ex &lh, const ex &rh)
81 {
82  construct_from_2_ex(lh,rh);
84 }
85 
87 {
90 }
91 
92 expairseq::expairseq(const epvector &v, const ex &oc, bool do_index_renaming)
93  : overall_coeff(oc)
94 {
95  GINAC_ASSERT(is_a<numeric>(oc));
96  construct_from_epvector(v, do_index_renaming);
97  GINAC_ASSERT(is_canonical());
98 }
99 
100 expairseq::expairseq(epvector && vp, const ex &oc, bool do_index_renaming)
101  : overall_coeff(oc)
102 {
103  GINAC_ASSERT(is_a<numeric>(oc));
104  construct_from_epvector(std::move(vp), do_index_renaming);
105  GINAC_ASSERT(is_canonical());
106 }
107 
109 // archiving
111 
112 void expairseq::read_archive(const archive_node &n, lst &sym_lst)
113 {
114  inherited::read_archive(n, sym_lst);
115  auto range = n.find_property_range("rest", "coeff");
116  seq.reserve((range.end-range.begin)/2);
117 
118  for (auto loc = range.begin; loc < range.end;) {
119  ex rest;
120  ex coeff;
121  n.find_ex_by_loc(loc++, rest, sym_lst);
122  n.find_ex_by_loc(loc++, coeff, sym_lst);
123  seq.emplace_back(expair(rest, coeff));
124  }
125 
126  n.find_ex("overall_coeff", overall_coeff, sym_lst);
127 
128  canonicalize();
129  GINAC_ASSERT(is_canonical());
130 }
131 
132 void expairseq::archive(archive_node &n) const
133 {
134  inherited::archive(n);
135  for (auto & i : seq) {
136  n.add_ex("rest", i.rest);
137  n.add_ex("coeff", i.coeff);
138  }
139  n.add_ex("overall_coeff", overall_coeff);
140 }
141 
142 
144 // functions overriding virtual functions from base classes
146 
147 // public
148 
149 void expairseq::do_print(const print_context & c, unsigned level) const
150 {
151  c.s << "[[";
152  printseq(c, ',', precedence(), level);
153  c.s << "]]";
154 }
155 
156 void expairseq::do_print_tree(const print_tree & c, unsigned level) const
157 {
158  c.s << std::string(level, ' ') << class_name() << " @" << this
159  << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
160  << ", nops=" << nops()
161  << std::endl;
162  size_t num = seq.size();
163  for (size_t i=0; i<num; ++i) {
164  seq[i].rest.print(c, level + c.delta_indent);
165  seq[i].coeff.print(c, level + c.delta_indent);
166  if (i != num - 1)
167  c.s << std::string(level + c.delta_indent, ' ') << "-----" << std::endl;
168  }
169  if (!overall_coeff.is_equal(default_overall_coeff())) {
170  c.s << std::string(level + c.delta_indent, ' ') << "-----" << std::endl
171  << std::string(level + c.delta_indent, ' ') << "overall_coeff" << std::endl;
172  overall_coeff.print(c, level + c.delta_indent);
173  }
174  c.s << std::string(level + c.delta_indent,' ') << "=====" << std::endl;
175 }
176 
177 bool expairseq::info(unsigned inf) const
178 {
179  switch(inf) {
180  case info_flags::expanded:
181  return (flags & status_flags::expanded);
182  case info_flags::has_indices: {
183  if (flags & status_flags::has_indices)
184  return true;
185  else if (flags & status_flags::has_no_indices)
186  return false;
187  for (auto & i : seq) {
188  if (i.rest.info(info_flags::has_indices)) {
189  this->setflag(status_flags::has_indices);
190  this->clearflag(status_flags::has_no_indices);
191  return true;
192  }
193  }
194  this->clearflag(status_flags::has_indices);
195  this->setflag(status_flags::has_no_indices);
196  return false;
197  }
198  }
199  return inherited::info(inf);
200 }
201 
202 size_t expairseq::nops() const
203 {
204  if (overall_coeff.is_equal(default_overall_coeff()))
205  return seq.size();
206  else
207  return seq.size()+1;
208 }
209 
210 ex expairseq::op(size_t i) const
211 {
212  if (i < seq.size())
213  return recombine_pair_to_ex(seq[i]);
214  GINAC_ASSERT(!overall_coeff.is_equal(default_overall_coeff()));
215  return overall_coeff;
216 }
217 
218 ex expairseq::map(map_function &f) const
219 {
220  epvector v;
221  v.reserve(seq.size()+1);
222 
223  for (auto & it : seq)
224  v.push_back(split_ex_to_pair(f(recombine_pair_to_ex(it))));
225 
226  if (overall_coeff.is_equal(default_overall_coeff()))
227  return thisexpairseq(std::move(v), default_overall_coeff(), true);
228  else {
229  ex newcoeff = f(overall_coeff);
230  if(is_a<numeric>(newcoeff))
231  return thisexpairseq(std::move(v), newcoeff, true);
232  else {
233  v.push_back(split_ex_to_pair(newcoeff));
234  return thisexpairseq(std::move(v), default_overall_coeff(), true);
235  }
236  }
237 }
238 
240 ex expairseq::eval() const
241 {
242  if (flags &status_flags::evaluated)
243  return *this;
244 
245  const epvector evaled = evalchildren();
246  if (!evaled.empty())
247  return dynallocate<expairseq>(std::move(evaled), overall_coeff).setflag(status_flags::evaluated);
248  else
249  return this->hold();
250 }
251 
253 {
254  epvector *newepv = nullptr;
255  for (auto i=epv.begin(); i!=epv.end(); ++i) {
256  if (newepv) {
257  newepv->push_back(i->conjugate());
258  continue;
259  }
260  expair x = i->conjugate();
261  if (x.is_equal(*i)) {
262  continue;
263  }
264  newepv = new epvector;
265  newepv->reserve(epv.size());
266  for (auto j=epv.begin(); j!=i; ++j) {
267  newepv->push_back(*j);
268  }
269  newepv->push_back(x);
270  }
271  return newepv;
272 }
273 
274 ex expairseq::conjugate() const
275 {
276  std::unique_ptr<epvector> newepv(conjugateepvector(seq));
277  ex x = overall_coeff.conjugate();
278  if (newepv) {
279  return thisexpairseq(std::move(*newepv), x);
280  }
281  if (are_ex_trivially_equal(x, overall_coeff)) {
282  return *this;
283  }
284  return thisexpairseq(seq, x);
285 }
286 
287 bool expairseq::match(const ex & pattern, exmap & repl_lst) const
288 {
289  // This differs from basic::match() because we want "a+b+c+d" to
290  // match "d+*+b" with "*" being "a+c", and we want to honor commutativity
291 
292  if (typeid(*this) == typeid(ex_to<basic>(pattern))) {
293 
294  // Check whether global wildcard (one that matches the "rest of the
295  // expression", like "*" above) is present
296  bool has_global_wildcard = false;
297  ex global_wildcard;
298  for (size_t i=0; i<pattern.nops(); i++) {
299  if (is_exactly_a<wildcard>(pattern.op(i))) {
300  has_global_wildcard = true;
301  global_wildcard = pattern.op(i);
302  break;
303  }
304  }
305 
306  // Even if the expression does not match the pattern, some of
307  // its subexpressions could match it. For example, x^5*y^(-1)
308  // does not match the pattern $0^5, but its subexpression x^5
309  // does. So, save repl_lst in order to not add bogus entries.
310  exmap tmp_repl = repl_lst;
311 
312  // Unfortunately, this is an O(N^2) operation because we can't
313  // sort the pattern in a useful way...
314 
315  // Chop into terms
316  exvector ops;
317  ops.reserve(nops());
318  for (size_t i=0; i<nops(); i++)
319  ops.push_back(op(i));
320 
321  // Now, for every term of the pattern, look for a matching term in
322  // the expression and remove the match
323  for (size_t i=0; i<pattern.nops(); i++) {
324  ex p = pattern.op(i);
325  if (has_global_wildcard && p.is_equal(global_wildcard))
326  continue;
327  auto it = ops.begin(), itend = ops.end();
328  while (it != itend) {
329  if (it->match(p, tmp_repl)) {
330  ops.erase(it);
331  goto found;
332  }
333  ++it;
334  }
335  return false; // no match found
336 found: ;
337  }
338 
339  if (has_global_wildcard) {
340 
341  // Assign all the remaining terms to the global wildcard (unless
342  // it has already been matched before, in which case the matches
343  // must be equal)
344  size_t num = ops.size();
345  epvector vp;
346  vp.reserve(num);
347  for (size_t i=0; i<num; i++)
348  vp.push_back(split_ex_to_pair(ops[i]));
349  ex rest = thisexpairseq(std::move(vp), default_overall_coeff());
350  for (auto & it : tmp_repl) {
351  if (it.first.is_equal(global_wildcard)) {
352  if (rest.is_equal(it.second)) {
353  repl_lst = tmp_repl;
354  return true;
355  }
356  return false;
357  }
358  }
359  repl_lst = tmp_repl;
360  repl_lst[global_wildcard] = rest;
361  return true;
362 
363  } else {
364 
365  // No global wildcard, then the match fails if there are any
366  // unmatched terms left
367  if (ops.empty()) {
368  repl_lst = tmp_repl;
369  return true;
370  }
371  return false;
372  }
373  }
374  return inherited::match(pattern, repl_lst);
375 }
376 
377 ex expairseq::subs(const exmap & m, unsigned options) const
378 {
379  epvector subsed = subschildren(m, options);
380  if (!subsed.empty())
381  return ex_to<basic>(thisexpairseq(std::move(subsed), overall_coeff, (options & subs_options::no_index_renaming) == 0));
382  else if ((options & subs_options::algebraic) && is_exactly_a<mul>(*this))
383  return static_cast<const mul *>(this)->algebraic_subs_mul(m, options);
384  else
385  return subs_one_level(m, options);
386 }
387 
388 // protected
389 
390 int expairseq::compare_same_type(const basic &other) const
391 {
392  GINAC_ASSERT(is_a<expairseq>(other));
393  const expairseq &o = static_cast<const expairseq &>(other);
394 
395  int cmpval;
396 
397  // compare number of elements
398  if (seq.size() != o.seq.size())
399  return (seq.size()<o.seq.size()) ? -1 : 1;
400 
401  // compare overall_coeff
402  cmpval = overall_coeff.compare(o.overall_coeff);
403  if (cmpval!=0)
404  return cmpval;
405 
406  auto cit1 = seq.begin(), last1 = seq.end();
407  auto cit2 = o.seq.begin(), last2 = o.seq.end();
408  for (; (cit1!=last1) && (cit2!=last2); ++cit1, ++cit2) {
409  cmpval = (*cit1).compare(*cit2);
410  if (cmpval!=0) return cmpval;
411  }
412 
413  GINAC_ASSERT(cit1==last1);
414  GINAC_ASSERT(cit2==last2);
415 
416  return 0;
417 }
418 
419 bool expairseq::is_equal_same_type(const basic &other) const
420 {
421  const expairseq &o = static_cast<const expairseq &>(other);
422 
423  // compare number of elements
424  if (seq.size()!=o.seq.size())
425  return false;
426 
427  // compare overall_coeff
428  if (!overall_coeff.is_equal(o.overall_coeff))
429  return false;
430 
431  auto cit2 = o.seq.begin();
432  for (auto & cit1 : seq) {
433  if (!cit1.is_equal(*cit2))
434  return false;
435  ++cit2;
436  }
437 
438  return true;
439 }
440 
441 unsigned expairseq::return_type() const
442 {
443  return return_types::noncommutative_composite;
444 }
445 
446 unsigned expairseq::calchash() const
447 {
448  unsigned v = make_hash_seed(typeid(*this));
449  for (auto & i : seq) {
450  v ^= i.rest.gethash();
451  v = rotate_left(v);
452  v ^= i.coeff.gethash();
453  }
454 
455  v ^= overall_coeff.gethash();
456 
457  // store calculated hash value only if object is already evaluated
458  if (flags &status_flags::evaluated) {
459  setflag(status_flags::hash_calculated);
460  hashvalue = v;
461  }
462 
463  return v;
464 }
465 
466 ex expairseq::expand(unsigned options) const
467 {
468  epvector expanded = expandchildren(options);
469  if (!expanded.empty()) {
470  return thisexpairseq(std::move(expanded), overall_coeff);
471  }
472  return (options == 0) ? setflag(status_flags::expanded) : *this;
473 }
474 
476 // new virtual functions which can be overridden by derived classes
478 
479 // protected
480 
489 ex expairseq::thisexpairseq(const epvector &v, const ex &oc, bool do_index_renaming) const
490 {
491  return expairseq(v, oc, do_index_renaming);
492 }
493 
494 ex expairseq::thisexpairseq(epvector && vp, const ex &oc, bool do_index_renaming) const
495 {
496  return expairseq(std::move(vp), oc, do_index_renaming);
497 }
498 
499 void expairseq::printpair(const print_context & c, const expair & p, unsigned upper_precedence) const
500 {
501  c.s << "[[";
502  p.rest.print(c, precedence());
503  c.s << ",";
504  p.coeff.print(c, precedence());
505  c.s << "]]";
506 }
507 
508 void expairseq::printseq(const print_context & c, char delim,
509  unsigned this_precedence,
510  unsigned upper_precedence) const
511 {
512  if (this_precedence <= upper_precedence)
513  c.s << "(";
514  auto it = seq.begin(), it_last = seq.end() - 1;
515  for (; it!=it_last; ++it) {
516  printpair(c, *it, this_precedence);
517  c.s << delim;
518  }
519  printpair(c, *it, this_precedence);
520  if (!overall_coeff.is_equal(default_overall_coeff())) {
521  c.s << delim;
522  overall_coeff.print(c, this_precedence);
523  }
524 
525  if (this_precedence <= upper_precedence)
526  c.s << ")";
527 }
528 
529 
532 expair expairseq::split_ex_to_pair(const ex &e) const
533 {
534  return expair(e,_ex1);
535 }
536 
537 
538 expair expairseq::combine_ex_with_coeff_to_pair(const ex &e,
539  const ex &c) const
540 {
541  GINAC_ASSERT(is_exactly_a<numeric>(c));
542 
543  return expair(e,c);
544 }
545 
546 
547 expair expairseq::combine_pair_with_coeff_to_pair(const expair &p,
548  const ex &c) const
549 {
550  GINAC_ASSERT(is_exactly_a<numeric>(p.coeff));
551  GINAC_ASSERT(is_exactly_a<numeric>(c));
552 
553  return expair(p.rest,ex_to<numeric>(p.coeff).mul_dyn(ex_to<numeric>(c)));
554 }
555 
556 
559 ex expairseq::recombine_pair_to_ex(const expair &p) const
560 {
561  return lst{p.rest, p.coeff};
562 }
563 
564 bool expairseq::expair_needs_further_processing(epp it)
565 {
566  if (is_exactly_a<numeric>(it->rest) &&
567  it->coeff.is_equal(_ex1)) {
568  // the pair {<n>, 1} has yet to be absorbed into overall_coeff
569  return true;
570  }
571  return false;
572 }
573 
574 ex expairseq::default_overall_coeff() const
575 {
576  return _ex0;
577 }
578 
579 void expairseq::combine_overall_coeff(const ex &c)
580 {
581  GINAC_ASSERT(is_exactly_a<numeric>(overall_coeff));
582  GINAC_ASSERT(is_exactly_a<numeric>(c));
583  overall_coeff = ex_to<numeric>(overall_coeff).add_dyn(ex_to<numeric>(c));
584 }
585 
586 void expairseq::combine_overall_coeff(const ex &c1, const ex &c2)
587 {
588  GINAC_ASSERT(is_exactly_a<numeric>(overall_coeff));
589  GINAC_ASSERT(is_exactly_a<numeric>(c1));
590  GINAC_ASSERT(is_exactly_a<numeric>(c2));
591  overall_coeff = ex_to<numeric>(overall_coeff).
592  add_dyn(ex_to<numeric>(c1).mul(ex_to<numeric>(c2)));
593 }
594 
595 bool expairseq::can_make_flat(const expair &p) const
596 {
597  return true;
598 }
599 
600 
602 // non-virtual functions in this class
604 
605 void expairseq::construct_from_2_ex(const ex &lh, const ex &rh)
606 {
607  const std::type_info& typeid_this = typeid(*this);
608  if (typeid(ex_to<basic>(lh)) == typeid_this) {
609  if (typeid(ex_to<basic>(rh)) == typeid_this) {
610  if (is_a<mul>(lh) && lh.info(info_flags::has_indices) &&
611  rh.info(info_flags::has_indices)) {
612  ex newrh=rename_dummy_indices_uniquely(lh, rh);
613  construct_from_2_expairseq(ex_to<expairseq>(lh),
614  ex_to<expairseq>(newrh));
615  }
616  else
617  construct_from_2_expairseq(ex_to<expairseq>(lh),
618  ex_to<expairseq>(rh));
619  return;
620  } else {
621  construct_from_expairseq_ex(ex_to<expairseq>(lh), rh);
622  return;
623  }
624  } else if (typeid(ex_to<basic>(rh)) == typeid_this) {
625  construct_from_expairseq_ex(ex_to<expairseq>(rh),lh);
626  return;
627  }
628 
629  if (is_exactly_a<numeric>(lh)) {
630  if (is_exactly_a<numeric>(rh)) {
631  combine_overall_coeff(lh);
632  combine_overall_coeff(rh);
633  } else {
634  combine_overall_coeff(lh);
635  seq.push_back(split_ex_to_pair(rh));
636  }
637  } else {
638  if (is_exactly_a<numeric>(rh)) {
639  combine_overall_coeff(rh);
640  seq.push_back(split_ex_to_pair(lh));
641  } else {
642  expair p1 = split_ex_to_pair(lh);
643  expair p2 = split_ex_to_pair(rh);
644 
645  int cmpval = p1.rest.compare(p2.rest);
646  if (cmpval==0) {
647  p1.coeff = ex_to<numeric>(p1.coeff).add_dyn(ex_to<numeric>(p2.coeff));
648  if (!ex_to<numeric>(p1.coeff).is_zero()) {
649  // no further processing is necessary, since this
650  // one element will usually be recombined in eval()
651  seq.push_back(p1);
652  }
653  } else {
654  seq.reserve(2);
655  if (cmpval<0) {
656  seq.push_back(p1);
657  seq.push_back(p2);
658  } else {
659  seq.push_back(p2);
660  seq.push_back(p1);
661  }
662  }
663  }
664  }
665 }
666 
667 void expairseq::construct_from_2_expairseq(const expairseq &s1,
668  const expairseq &s2)
669 {
670  combine_overall_coeff(s1.overall_coeff);
671  combine_overall_coeff(s2.overall_coeff);
672 
673  auto first1 = s1.seq.begin(), last1 = s1.seq.end();
674  auto first2 = s2.seq.begin(), last2 = s2.seq.end();
675 
676  seq.reserve(s1.seq.size()+s2.seq.size());
677 
678  bool needs_further_processing=false;
679 
680  while (first1!=last1 && first2!=last2) {
681  int cmpval = (*first1).rest.compare((*first2).rest);
682 
683  if (cmpval==0) {
684  // combine terms
685  const numeric &newcoeff = ex_to<numeric>(first1->coeff).
686  add(ex_to<numeric>(first2->coeff));
687  if (!newcoeff.is_zero()) {
688  seq.push_back(expair(first1->rest,newcoeff));
689  if (expair_needs_further_processing(seq.end()-1)) {
690  needs_further_processing = true;
691  }
692  }
693  ++first1;
694  ++first2;
695  } else if (cmpval<0) {
696  seq.push_back(*first1);
697  ++first1;
698  } else {
699  seq.push_back(*first2);
700  ++first2;
701  }
702  }
703 
704  while (first1!=last1) {
705  seq.push_back(*first1);
706  ++first1;
707  }
708  while (first2!=last2) {
709  seq.push_back(*first2);
710  ++first2;
711  }
712 
713  if (needs_further_processing) {
714  // Clear seq and start over.
715  epvector v = std::move(seq);
716  construct_from_epvector(std::move(v));
717  }
718 }
719 
720 void expairseq::construct_from_expairseq_ex(const expairseq &s,
721  const ex &e)
722 {
723  combine_overall_coeff(s.overall_coeff);
724  if (is_exactly_a<numeric>(e)) {
725  combine_overall_coeff(e);
726  seq = s.seq;
727  return;
728  }
729 
730  auto first = s.seq.begin(), last = s.seq.end();
731  expair p = split_ex_to_pair(e);
732 
733  seq.reserve(s.seq.size()+1);
734  bool p_pushed = false;
735 
736  bool needs_further_processing=false;
737 
738  // merge p into s.seq
739  while (first!=last) {
740  int cmpval = (*first).rest.compare(p.rest);
741  if (cmpval==0) {
742  // combine terms
743  const numeric &newcoeff = ex_to<numeric>(first->coeff).
744  add(ex_to<numeric>(p.coeff));
745  if (!newcoeff.is_zero()) {
746  seq.push_back(expair(first->rest,newcoeff));
747  if (expair_needs_further_processing(seq.end()-1))
748  needs_further_processing = true;
749  }
750  ++first;
751  p_pushed = true;
752  break;
753  } else if (cmpval<0) {
754  seq.push_back(*first);
755  ++first;
756  } else {
757  seq.push_back(p);
758  p_pushed = true;
759  break;
760  }
761  }
762 
763  if (p_pushed) {
764  // while loop exited because p was pushed, now push rest of s.seq
765  while (first!=last) {
766  seq.push_back(*first);
767  ++first;
768  }
769  } else {
770  // while loop exited because s.seq was pushed, now push p
771  seq.push_back(p);
772  }
773 
774  if (needs_further_processing) {
775  // Clear seq and start over.
776  epvector v = std::move(seq);
777  construct_from_epvector(std::move(v));
778  }
779 }
780 
781 void expairseq::construct_from_exvector(const exvector &v)
782 {
783  // simplifications: +(a,+(b,c),d) -> +(a,b,c,d) (associativity)
784  // +(d,b,c,a) -> +(a,b,c,d) (canonicalization)
785  // +(...,x,*(x,c1),*(x,c2)) -> +(...,*(x,1+c1+c2)) (c1, c2 numeric)
786  // (same for (+,*) -> (*,^)
787 
788  make_flat(v);
789  canonicalize();
790  combine_same_terms_sorted_seq();
791 }
792 
793 void expairseq::construct_from_epvector(const epvector &v, bool do_index_renaming)
794 {
795  // simplifications: +(a,+(b,c),d) -> +(a,b,c,d) (associativity)
796  // +(d,b,c,a) -> +(a,b,c,d) (canonicalization)
797  // +(...,x,*(x,c1),*(x,c2)) -> +(...,*(x,1+c1+c2)) (c1, c2 numeric)
798  // same for (+,*) -> (*,^)
799 
800  make_flat(v, do_index_renaming);
801  canonicalize();
802  combine_same_terms_sorted_seq();
803 }
804 
805 void expairseq::construct_from_epvector(epvector &&v, bool do_index_renaming)
806 {
807  // simplifications: +(a,+(b,c),d) -> +(a,b,c,d) (associativity)
808  // +(d,b,c,a) -> +(a,b,c,d) (canonicalization)
809  // +(...,x,*(x,c1),*(x,c2)) -> +(...,*(x,1+c1+c2)) (c1, c2 numeric)
810  // same for (+,*) -> (*,^)
811 
812  make_flat(std::move(v), do_index_renaming);
813  canonicalize();
814  combine_same_terms_sorted_seq();
815 }
816 
819 void expairseq::make_flat(const exvector &v)
820 {
821  // count number of operands which are of same expairseq derived type
822  // and their cumulative number of operands
823  int nexpairseqs = 0;
824  int noperands = 0;
825  bool do_idx_rename = false;
826 
827  const std::type_info& typeid_this = typeid(*this);
828  for (auto & cit : v) {
829  if (typeid(ex_to<basic>(cit)) == typeid_this) {
830  ++nexpairseqs;
831  noperands += ex_to<expairseq>(cit).seq.size();
832  }
833  if (is_a<mul>(*this) && (!do_idx_rename) &&
834  cit.info(info_flags::has_indices))
835  do_idx_rename = true;
836  }
837 
838  // reserve seq and coeffseq which will hold all operands
839  seq.reserve(v.size()+noperands-nexpairseqs);
840 
841  // copy elements and split off numerical part
842  make_flat_inserter mf(v, do_idx_rename);
843  for (auto & cit : v) {
844  if (typeid(ex_to<basic>(cit)) == typeid_this) {
845  ex newfactor = mf.handle_factor(cit, _ex1);
846  const expairseq &subseqref = ex_to<expairseq>(newfactor);
847  combine_overall_coeff(subseqref.overall_coeff);
848  for (auto & cit_s : subseqref.seq) {
849  seq.push_back(cit_s);
850  }
851  } else {
852  if (is_exactly_a<numeric>(cit))
853  combine_overall_coeff(cit);
854  else {
855  ex newfactor = mf.handle_factor(cit, _ex1);
856  seq.push_back(split_ex_to_pair(newfactor));
857  }
858  }
859  }
860 }
861 
864 void expairseq::make_flat(const epvector &v, bool do_index_renaming)
865 {
866  // count number of operands which are of same expairseq derived type
867  // and their cumulative number of operands
868  int nexpairseqs = 0;
869  int noperands = 0;
870  bool really_need_rename_inds = false;
871 
872  const std::type_info& typeid_this = typeid(*this);
873  for (auto & cit : v) {
874  if (typeid(ex_to<basic>(cit.rest)) == typeid_this) {
875  ++nexpairseqs;
876  noperands += ex_to<expairseq>(cit.rest).seq.size();
877  }
878  if ((!really_need_rename_inds) && is_a<mul>(*this) &&
879  cit.rest.info(info_flags::has_indices))
880  really_need_rename_inds = true;
881  }
882  do_index_renaming = do_index_renaming && really_need_rename_inds;
883 
884  // reserve seq and coeffseq which will hold all operands
885  seq.reserve(v.size()+noperands-nexpairseqs);
886  make_flat_inserter mf(v, do_index_renaming);
887 
888  // copy elements and split off numerical part
889  for (auto & cit : v) {
890  if (typeid(ex_to<basic>(cit.rest)) == typeid_this &&
891  this->can_make_flat(cit)) {
892  ex newrest = mf.handle_factor(cit.rest, cit.coeff);
893  const expairseq &subseqref = ex_to<expairseq>(newrest);
894  combine_overall_coeff(subseqref.overall_coeff, cit.coeff);
895  for (auto & cit_s : subseqref.seq) {
896  seq.push_back(expair(cit_s.rest,
897  ex_to<numeric>(cit_s.coeff).mul_dyn(ex_to<numeric>(cit.coeff))));
898  }
899  } else {
900  if (cit.is_canonical_numeric())
901  combine_overall_coeff(mf.handle_factor(cit.rest, _ex1));
902  else {
903  ex rest = cit.rest;
904  ex newrest = mf.handle_factor(rest, cit.coeff);
905  if (are_ex_trivially_equal(newrest, rest))
906  seq.push_back(cit);
907  else
908  seq.push_back(expair(newrest, cit.coeff));
909  }
910  }
911  }
912 }
913 
916 {
917  std::sort(seq.begin(), seq.end(), expair_rest_is_less());
918 }
919 
920 
924 void expairseq::combine_same_terms_sorted_seq()
925 {
926  if (seq.size()<2)
927  return;
928 
929  bool needs_further_processing = false;
930 
931  auto itin1 = seq.begin();
932  auto itin2 = itin1 + 1;
933  auto itout = itin1;
934  auto last = seq.end();
935  // must_copy will be set to true the first time some combination is
936  // possible from then on the sequence has changed and must be compacted
937  bool must_copy = false;
938  while (itin2!=last) {
939  if (itin1->rest.compare(itin2->rest)==0) {
940  itin1->coeff = ex_to<numeric>(itin1->coeff).
941  add_dyn(ex_to<numeric>(itin2->coeff));
942  if (expair_needs_further_processing(itin1))
943  needs_further_processing = true;
944  must_copy = true;
945  } else {
946  if (!ex_to<numeric>(itin1->coeff).is_zero()) {
947  if (must_copy)
948  *itout = *itin1;
949  ++itout;
950  }
951  itin1 = itin2;
952  }
953  ++itin2;
954  }
955  if (!ex_to<numeric>(itin1->coeff).is_zero()) {
956  if (must_copy)
957  *itout = *itin1;
958  ++itout;
959  }
960  if (itout!=last)
961  seq.erase(itout,last);
962 
963  if (needs_further_processing) {
964  // Clear seq and start over.
965  epvector v = std::move(seq);
966  construct_from_epvector(std::move(v));
967  }
968 }
969 
972 bool expairseq::is_canonical() const
973 {
974  if (seq.size() <= 1)
975  return 1;
976 
977  auto it = seq.begin(), itend = seq.end();
978  auto it_last = it;
979  for (++it; it!=itend; it_last=it, ++it) {
980  if (!(it_last->is_less(*it) || it_last->is_equal(*it))) {
981  if (!is_exactly_a<numeric>(it_last->rest) ||
982  !is_exactly_a<numeric>(it->rest)) {
983  // double test makes it easier to set a breakpoint...
984  if (!is_exactly_a<numeric>(it_last->rest) ||
985  !is_exactly_a<numeric>(it->rest)) {
986  printpair(std::clog, *it_last, 0);
987  std::clog << ">";
988  printpair(std::clog, *it, 0);
989  std::clog << "\n";
990  std::clog << "pair1:" << std::endl;
991  it_last->rest.print(print_tree(std::clog));
992  it_last->coeff.print(print_tree(std::clog));
993  std::clog << "pair2:" << std::endl;
994  it->rest.print(print_tree(std::clog));
995  it->coeff.print(print_tree(std::clog));
996  return 0;
997  }
998  }
999  }
1000  }
1001  return 1;
1002 }
1003 
1009 epvector expairseq::expandchildren(unsigned options) const
1010 {
1011  auto cit = seq.begin(), last = seq.end();
1012  while (cit!=last) {
1013  const ex expanded_ex = cit->rest.expand(options);
1014  if (!are_ex_trivially_equal(cit->rest,expanded_ex)) {
1015 
1016  // something changed, copy seq, eval and return it
1017  epvector s;
1018  s.reserve(seq.size());
1019 
1020  // copy parts of seq which are known not to have changed
1021  s.insert(s.begin(), seq.begin(), cit);
1022 
1023  // copy first changed element
1024  s.push_back(expair(expanded_ex, cit->coeff));
1025  ++cit;
1026 
1027  // copy rest
1028  while (cit != last) {
1029  s.push_back(expair(cit->rest.expand(options), cit->coeff));
1030  ++cit;
1031  }
1032  return s;
1033  }
1034 
1035  ++cit;
1036  }
1037 
1038  return epvector(); // empty signalling nothing has changed
1039 }
1040 
1041 
1047 epvector expairseq::evalchildren() const
1048 {
1049  auto cit = seq.begin(), last = seq.end();
1050  while (cit!=last) {
1051  const expair evaled_pair = combine_ex_with_coeff_to_pair(cit->rest, cit->coeff);
1052  if (unlikely(!evaled_pair.is_equal(*cit))) {
1053 
1054  // something changed: copy seq, eval and return it
1055  epvector s;
1056  s.reserve(seq.size());
1057 
1058  // copy parts of seq which are known not to have changed
1059  s.insert(s.begin(), seq.begin(), cit);
1060 
1061  // copy first changed element
1062  s.push_back(evaled_pair);
1063  ++cit;
1064 
1065  // copy rest
1066  while (cit != last) {
1067  s.push_back(combine_ex_with_coeff_to_pair(cit->rest, cit->coeff));
1068  ++cit;
1069  }
1070  return s;
1071  }
1072 
1073  ++cit;
1074  }
1075 
1076  return epvector(); // signalling nothing has changed
1077 }
1078 
1084 epvector expairseq::subschildren(const exmap & m, unsigned options) const
1085 {
1086  // When any of the objects to be substituted is a product or power
1087  // we have to recombine the pairs because the numeric coefficients may
1088  // be part of the search pattern.
1089  if (!(options & (subs_options::pattern_is_product | subs_options::pattern_is_not_product))) {
1090 
1091  // Search the list of substitutions and cache our findings
1092  for (auto & it : m) {
1093  if (is_exactly_a<mul>(it.first) || is_exactly_a<power>(it.first)) {
1094  options |= subs_options::pattern_is_product;
1095  break;
1096  }
1097  }
1098  if (!(options & subs_options::pattern_is_product))
1099  options |= subs_options::pattern_is_not_product;
1100  }
1101 
1102  if (options & subs_options::pattern_is_product) {
1103 
1104  // Substitute in the recombined pairs
1105  auto cit = seq.begin(), last = seq.end();
1106  while (cit != last) {
1107 
1108  const ex &orig_ex = recombine_pair_to_ex(*cit);
1109  const ex &subsed_ex = orig_ex.subs(m, options);
1110  if (!are_ex_trivially_equal(orig_ex, subsed_ex)) {
1111 
1112  // Something changed: copy seq, subs and return it
1113  epvector s;
1114  s.reserve(seq.size());
1115 
1116  // Copy parts of seq which are known not to have changed
1117  s.insert(s.begin(), seq.begin(), cit);
1118 
1119  // Copy first changed element
1120  s.push_back(split_ex_to_pair(subsed_ex));
1121  ++cit;
1122 
1123  // Copy rest
1124  while (cit != last) {
1125  s.push_back(split_ex_to_pair(recombine_pair_to_ex(*cit).subs(m, options)));
1126  ++cit;
1127  }
1128  return s;
1129  }
1130 
1131  ++cit;
1132  }
1133 
1134  } else {
1135 
1136  // Substitute only in the "rest" part of the pairs
1137  auto cit = seq.begin(), last = seq.end();
1138  while (cit != last) {
1139 
1140  const ex subsed_rest = cit->rest.subs(m, options);
1141  const expair subsed_pair = combine_ex_with_coeff_to_pair(subsed_rest, cit->coeff);
1142  if (!subsed_pair.is_equal(*cit)) {
1143 
1144  // Something changed: copy seq, subs and return it
1145  epvector s;
1146  s.reserve(seq.size());
1147 
1148  // Copy parts of seq which are known not to have changed
1149  s.insert(s.begin(), seq.begin(), cit);
1150 
1151  // Copy first changed element
1152  s.push_back(subsed_pair);
1153  ++cit;
1154 
1155  // Copy rest
1156  while (cit != last) {
1157  s.push_back(combine_ex_with_coeff_to_pair(cit->rest.subs(m, options), cit->coeff));
1158  ++cit;
1159  }
1160  return s;
1161  }
1162 
1163  ++cit;
1164  }
1165  }
1166 
1167  // Nothing has changed
1168  return epvector();
1169 }
1170 
1172 // static member variables
1174 
1175 } // namespace GiNaC
size_t nops(const ex &thisex)
Definition: ex.h:712
Interface to GiNaC&#39;s symbolic exponentiation (basis^exponent).
void do_print(const print_context &c, unsigned level) const
Interface to GiNaC&#39;s indexed expressions.
int canonicalize(exvector::iterator v, const symmetry &symm)
Canonicalize the order of elements of an expression vector, according to the symmetry properties defi...
Definition: symmetry.cpp:438
Definition: add.cpp:38
Archiving of GiNaC expressions.
Interface to GiNaC&#39;s sums of expressions.
Interface to sequences of expression pairs.
void do_print_tree(const print_tree &c, unsigned level) const
const ex _ex0
Definition: utils.cpp:177
bool are_ex_trivially_equal(const ex &e1, const ex &e2)
Compare two objects of class quickly without doing a deep tree traversal.
Definition: ex.h:684
ex x
Definition: factor.cpp:1641
Interface to GiNaC&#39;s products of expressions.
Interface to GiNaC&#39;s wildcard objects.
container< std::list > lst
Definition: lst.h:32
Interface to several small and furry utilities needed within GiNaC but not of any interest to the use...
const ex _ex1
Definition: utils.cpp:193
Definition of optimizing macros.
#define unlikely(cond)
Definition: compiler.h:31
ex conjugate(const ex &thisex)
Definition: ex.h:718
print_func< print_context >(&varidx::do_print). print_func< print_latex >(&varidx
Definition: idx.cpp:45
unsigned options
Definition: factor.cpp:2480
mvec m
Definition: factor.cpp:771
#define GINAC_ASSERT(X)
Assertion macro for checking invariances.
Definition: assertion.h:33
bool match(const ex &thisex, const ex &pattern, exmap &repl_lst)
Definition: ex.h:784
static unsigned make_hash_seed(const std::type_info &tinfo)
We need a hash function which gives different values for objects of different types.
Definition: hash_seed.h:36
epvector * conjugateepvector(const epvector &)
Complex conjugate every element of an epvector.
void construct_from_2_ex(const ex &lh, const ex &rh)
GINAC_IMPLEMENT_REGISTERED_CLASS_OPT(add, expairseq, print_func< print_context >(&add::do_print). print_func< print_latex >(&add::do_print_latex). print_func< print_csrc >(&add::do_print_csrc). print_func< print_tree >(&add::do_print_tree). print_func< print_python_repr >(&add::do_print_python_repr)) add
Definition: add.cpp:40
ex op(const ex &thisex, size_t i)
Definition: ex.h:811
std::vector< ex > exvector
Definition: basic.h:46
std::map< ex, ex, ex_is_less > exmap
Definition: basic.h:50
epvector::iterator epp
expair-vector pointer
Definition: expairseq.h:34
Interface to GiNaC&#39;s overloaded operators.
size_t last
Definition: factor.cpp:1465
size_t n
Definition: factor.cpp:1463
Definition of GiNaC&#39;s lst.
unsigned rotate_left(unsigned n)
Rotate bits of unsigned value by one bit to the left.
Definition: utils.h:48
lst rename_dummy_indices_uniquely(const exvector &va, const exvector &vb)
Similar to above, where va and vb are the same and the return value is a list of two lists for substi...
Definition: indexed.cpp:1460
ex coeff(const ex &thisex, const ex &s, int n=1)
Definition: ex.h:742
expairseq(const ex &lh, const ex &rh)
Interface to relations between expressions.
std::vector< expair > epvector
expair-vector
Definition: expairseq.h:33
bool is_canonical() const
size_t c
Definition: factor.cpp:770
void construct_from_exvector(const exvector &v)
ex eval(const ex &thisex)
Definition: ex.h:766
ex subs(const ex &thisex, const exmap &m, unsigned options=0)
Definition: ex.h:831
ex expand(const ex &thisex, unsigned options=0)
Definition: ex.h:715
Type-specific hash seed.

This page is part of the GiNaC developer's reference. It was generated automatically by doxygen. For an introduction, see the tutorial.