Needed to call get_all_dummy_indices_safely more often.

[ginac.git] / ginac / mul.cpp

diff --git a/ginac/mul.cpp b/ginac/mul.cpp
index 962947898baf696c245ab816866fb286a39c5fdf..4cea8167b4ee419081c67f275678c1e78b3565df 100644 (file)
--- a/ginac/mul.cpp
+++ b/ginac/mul.cpp
@@ -27,8 +27,6 @@
 
 #include "mul.h"
 #include "add.h"
-#include "color.h"
-#include "clifford.h"
 #include "power.h"
 #include "operators.h"
 #include "matrix.h"
@@ -463,6 +461,41 @@ ex mul::evalf(int level) const
        return mul(s, overall_coeff.evalf(level));
 }
 
+void mul::find_real_imag(ex & rp, ex & ip) const
+{
+       rp = overall_coeff.real_part();
+       ip = overall_coeff.imag_part();
+       for (epvector::const_iterator i=seq.begin(); i!=seq.end(); ++i) {
+               ex factor = recombine_pair_to_ex(*i);
+               ex new_rp = factor.real_part();
+               ex new_ip = factor.imag_part();
+               if(new_ip.is_zero()) {
+                       rp *= new_rp;
+                       ip *= new_rp;
+               } else {
+                       ex temp = rp*new_rp - ip*new_ip;
+                       ip = ip*new_rp + rp*new_ip;
+                       rp = temp;
+               }
+       }
+       rp = rp.expand();
+       ip = ip.expand();
+}
+
+ex mul::real_part() const
+{
+       ex rp, ip;
+       find_real_imag(rp, ip);
+       return rp;
+}
+
+ex mul::imag_part() const
+{
+       ex rp, ip;
+       find_real_imag(rp, ip);
+       return ip;
+}
+
 ex mul::evalm() const
 {
        // numeric*matrix
@@ -597,6 +630,22 @@ bool algebraic_match_mul_with_mul(const mul &e, const ex &pat, lst &repls,
        return false;
 }
 
+bool mul::has(const ex & pattern, unsigned options) const
+{
+       if(!(options&has_options::algebraic))
+               return basic::has(pattern,options);
+       if(is_a<mul>(pattern)) {
+               lst repls;
+               int nummatches = std::numeric_limits<int>::max();
+               std::vector<bool> subsed(seq.size(), false);
+               std::vector<bool> matched(seq.size(), false);
+               if(algebraic_match_mul_with_mul(*this, pattern, repls, 0, nummatches,
+                               subsed, matched))
+                       return true;
+       }
+       return basic::has(pattern, options);
+}
+
 ex mul::algebraic_subs_mul(const exmap & m, unsigned options) const
 {      
        std::vector<bool> subsed(seq.size(), false);
@@ -716,21 +765,8 @@ unsigned mul::return_type() const
                }
                if ((rt == return_types::noncommutative) && (!all_commutative)) {
                        // another nc element found, compare type_infos
-                       if (noncommutative_element->rest.return_type_tinfo()->tinfo() == &clifford::tinfo_static) {
-                               if (i->rest.return_type_tinfo()->tinfo() != &clifford::tinfo_static ||
-                                   ((clifford*)(noncommutative_element->rest.return_type_tinfo()))->get_representation_label() !=
-                                   ((clifford*)(i->rest.return_type_tinfo()))->get_representation_label()) {
-                                       // diffent types -> mul is ncc
-                                       return return_types::noncommutative_composite;
-                               }
-                       } else if (noncommutative_element->rest.return_type_tinfo()->tinfo() == &color::tinfo_static) {
-                               if (i->rest.return_type_tinfo()->tinfo() != &color::tinfo_static ||
-                                   ((color*)(noncommutative_element->rest.return_type_tinfo()))->get_representation_label() !=
-                                   ((color*)(i->rest.return_type_tinfo()))->get_representation_label()) {
-                                       // diffent types -> mul is ncc
-                                       return return_types::noncommutative_composite;
-                               }
-                       } else if (noncommutative_element->rest.return_type_tinfo()->tinfo() != i->rest.return_type_tinfo()->tinfo()) {
+                       if (noncommutative_element->rest.return_type_tinfo() != i->rest.return_type_tinfo()) {
+                                       // different types -> mul is ncc
                                        return return_types::noncommutative_composite;
                        }
                }
@@ -740,7 +776,7 @@ unsigned mul::return_type() const
        return all_commutative ? return_types::commutative : return_types::noncommutative;
 }
    
-const basic* mul::return_type_tinfo() const
+tinfo_t mul::return_type_tinfo() const
 {
        if (seq.empty())
                return this;  // mul without factors: should not happen
@@ -935,11 +971,11 @@ ex mul::expand(unsigned options) const
                                exvector add1_dummy_indices, add2_dummy_indices, add_indices;
 
                                for (epvector::const_iterator i=add1begin; i!=add1end; ++i) {
-                                       add_indices = get_all_dummy_indices(i->rest);
+                                       add_indices = get_all_dummy_indices_safely(i->rest);
                                        add1_dummy_indices.insert(add1_dummy_indices.end(), add_indices.begin(), add_indices.end());
                                }
                                for (epvector::const_iterator i=add2begin; i!=add2end; ++i) {
-                                       add_indices = get_all_dummy_indices(i->rest);
+                                       add_indices = get_all_dummy_indices_safely(i->rest);
                                        add2_dummy_indices.insert(add2_dummy_indices.end(), add_indices.begin(), add_indices.end());
                                }
 
@@ -986,7 +1022,7 @@ ex mul::expand(unsigned options) const
                size_t n = last_expanded.nops();
                exvector distrseq;
                distrseq.reserve(n);
-               exvector va = get_all_dummy_indices(mul(non_adds));
+               exvector va = get_all_dummy_indices_safely(mul(non_adds));
                sort(va.begin(), va.end(), ex_is_less());
 
                for (size_t i=0; i<n; ++i) {