- collect_common_factors() works better with negative powers

- simplify_indexed() converts "gamma~mu*p.mu" to "p\"
- fixed a bug in canonicalize_clifford() when gammaL/R were in a string

diff --git a/ginac/clifford.cpp b/ginac/clifford.cpp
index 0a10d4525854f5f575297234aefa33ba0f3fe2db..af8b749c34a453fadff041722785ca71fd0cd1ab 100644 (file)
--- a/ginac/clifford.cpp
+++ b/ginac/clifford.cpp
@@ -208,6 +208,7 @@ bool diracgamma::contract_with(exvector::iterator self, exvector::iterator other
        GINAC_ASSERT(is_a<indexed>(*other));
        GINAC_ASSERT(is_a<diracgamma>(self->op(0)));
        unsigned char rl = ex_to<clifford>(*self).get_representation_label();
+       ex dim = ex_to<idx>(self->op(1)).get_dim();
 
        if (is_a<clifford>(*other)) {
 
@@ -215,8 +216,6 @@ bool diracgamma::contract_with(exvector::iterator self, exvector::iterator other
                if (ex_to<clifford>(*other).get_representation_label() != rl)
                        return false;
 
-               ex dim = ex_to<idx>(self->op(1)).get_dim();
-
                // gamma~mu gamma.mu = dim ONE
                if (other - self == 1) {
                        *self = dim;
@@ -278,6 +277,13 @@ bool diracgamma::contract_with(exvector::iterator self, exvector::iterator other
                        *other = _ex1;
                        return true;
                }
+
+       } else if (is_a<symbol>(other->op(0)) && other->nops() == 2) {
+
+               // x.mu gamma~mu -> x-slash
+               *self = dirac_slash(other->op(0), dim, rl);
+               *other = _ex1;
+               return true;
        }
 
        return false;
@@ -718,7 +724,7 @@ ex canonicalize_clifford(const ex & e)
 
                        // Stupid recursive bubble sort because we only want to swap adjacent gammas
                        exvector::iterator it = v.begin(), next_to_last = v.end() - 1;
-                       if (is_a<diracgamma5>(it->op(0)))
+                       if (is_a<diracgamma5>(it->op(0)) || is_a<diracgammaL>(it->op(0)) || is_a<diracgammaR>(it->op(0)))
                                ++it;
                        while (it != next_to_last) {
                                if (it[0].compare(it[1]) > 0) {

diff --git a/ginac/normal.cpp b/ginac/normal.cpp
index 67091f5003aef503359fdee4814f3c4e07d504be..e1d7349fa509b600cf8321c54a5ca1ed0a93180d 100644 (file)
--- a/ginac/normal.cpp
+++ b/ginac/normal.cpp
@@ -1741,6 +1741,7 @@ static exvector sqrfree_yun(const ex &a, const symbol &x)
        return res;
 }
 
+
 /** Compute a square-free factorization of a multivariate polynomial in Q[X].
  *
  *  @param a  multivariate polynomial over Q[X]
@@ -1842,6 +1843,7 @@ ex sqrfree(const ex &a, const lst &l)
        return result * lcm.inverse();
 }
 
+
 /** Compute square-free partial fraction decomposition of rational function
  *  a(x).
  *
@@ -1912,96 +1914,6 @@ ex sqrfree_parfrac(const ex & a, const symbol & x)
 }
 
 
-/** Remove the common factor in the terms of a sum 'e' by calculating the GCD,
- *  and multiply it into the expression 'factor' (which needs to be initialized
- *  to 1, unless you're accumulating factors). */
-static ex find_common_factor(const ex & e, ex & factor)
-{
-       if (is_a<add>(e)) {
-
-               unsigned num = e.nops();
-               exvector terms; terms.reserve(num);
-               lst repl;
-               ex gc;
-
-               // Find the common GCD
-               for (unsigned i=0; i<num; i++) {
-                       ex x = e.op(i).to_rational(repl);
-                       if (is_a<add>(x) || is_a<mul>(x)) {
-                               ex f = 1;
-                               x = find_common_factor(x, f);
-                               x *= f;
-                       }
-
-                       if (i == 0)
-                               gc = x;
-                       else
-                               gc = gcd(gc, x);
-
-                       terms.push_back(x);
-               }
-
-               if (gc.is_equal(_ex1))
-                       return e;
-
-               // The GCD is the factor we pull out
-               factor *= gc.subs(repl);
-
-               // Now divide all terms by the GCD
-               for (unsigned i=0; i<num; i++) {
-                       ex x;
-
-                       // Try to avoid divide() because it expands the polynomial
-                       ex &t = terms[i];
-                       if (is_a<mul>(t)) {
-                               for (unsigned j=0; j<t.nops(); j++) {
-                                       if (t.op(j).is_equal(gc)) {
-                                               exvector v; v.reserve(t.nops());
-                                               for (unsigned k=0; k<t.nops(); k++) {
-                                                       if (k == j)
-                                                               v.push_back(_ex1);
-                                                       else
-                                                               v.push_back(t.op(k));
-                                               }
-                                               t = (new mul(v))->setflag(status_flags::dynallocated).subs(repl);
-                                               goto term_done;
-                                       }
-                               }
-                       }
-
-                       divide(t, gc, x);
-                       t = x.subs(repl);
-term_done:     ;
-               }
-               return (new add(terms))->setflag(status_flags::dynallocated);
-
-       } else if (is_a<mul>(e)) {
-
-               exvector v;
-               for (unsigned i=0; i<e.nops(); i++)
-                       v.push_back(find_common_factor(e.op(i), factor));
-               return (new mul(v))->setflag(status_flags::dynallocated);
-
-       } else
-               return e;
-}
-
-
-/** Collect common factors in sums. This converts expressions like
- *  'a*(b*x+b*y)' to 'a*b*(x+y)'. */
-ex collect_common_factors(const ex & e)
-{
-       if (is_a<add>(e) || is_a<mul>(e)) {
-
-               ex factor = 1;
-               ex r = find_common_factor(e, factor);
-               return factor * r;
-
-       } else
-               return e;
-}
-
-
 /*
  *  Normal form of rational functions
  */
@@ -2503,4 +2415,106 @@ ex expairseq::to_rational(lst &repl_lst) const
 }
 
 
+/** Remove the common factor in the terms of a sum 'e' by calculating the GCD,
+ *  and multiply it into the expression 'factor' (which needs to be initialized
+ *  to 1, unless you're accumulating factors). */
+static ex find_common_factor(const ex & e, ex & factor, lst & repl)
+{
+       if (is_a<add>(e)) {
+
+               unsigned num = e.nops();
+               exvector terms; terms.reserve(num);
+               ex gc;
+
+               // Find the common GCD
+               for (unsigned i=0; i<num; i++) {
+                       ex x = e.op(i).to_rational(repl);
+
+                       if (is_a<add>(x) || is_a<mul>(x)) {
+                               ex f = 1;
+                               x = find_common_factor(x, f, repl);
+                               x *= f;
+                       }
+
+                       if (i == 0)
+                               gc = x;
+                       else
+                               gc = gcd(gc, x);
+
+                       terms.push_back(x);
+               }
+
+               if (gc.is_equal(_ex1))
+                       return e;
+
+               // The GCD is the factor we pull out
+               factor *= gc;
+
+               // Now divide all terms by the GCD
+               for (unsigned i=0; i<num; i++) {
+                       ex x;
+
+                       // Try to avoid divide() because it expands the polynomial
+                       ex &t = terms[i];
+                       if (is_a<mul>(t)) {
+                               for (unsigned j=0; j<t.nops(); j++) {
+                                       if (t.op(j).is_equal(gc)) {
+                                               exvector v; v.reserve(t.nops());
+                                               for (unsigned k=0; k<t.nops(); k++) {
+                                                       if (k == j)
+                                                               v.push_back(_ex1);
+                                                       else
+                                                               v.push_back(t.op(k));
+                                               }
+                                               t = (new mul(v))->setflag(status_flags::dynallocated);
+                                               goto term_done;
+                                       }
+                               }
+                       }
+
+                       divide(t, gc, x);
+                       t = x;
+term_done:     ;
+               }
+               return (new add(terms))->setflag(status_flags::dynallocated);
+
+       } else if (is_a<mul>(e)) {
+
+               unsigned num = e.nops();
+               exvector v; v.reserve(num);
+
+               for (unsigned i=0; i<num; i++)
+                       v.push_back(find_common_factor(e.op(i), factor, repl));
+
+               return (new mul(v))->setflag(status_flags::dynallocated);
+
+       } else if (is_a<power>(e)) {
+
+               ex x = e.to_rational(repl);
+               if (is_a<power>(x) && x.op(1).info(info_flags::negative))
+                       return replace_with_symbol(x, repl);
+               else
+                       return x;
+
+       } else
+               return e;
+}
+
+
+/** Collect common factors in sums. This converts expressions like
+ *  'a*(b*x+b*y)' to 'a*b*(x+y)'. */
+ex collect_common_factors(const ex & e)
+{
+       if (is_a<add>(e) || is_a<mul>(e)) {
+
+               lst repl;
+               ex factor = 1;
+               ex r = find_common_factor(e, factor, repl);
+               return factor.subs(repl) * r.subs(repl);
+
+       } else
+               return e;
+}
+
+
 } // namespace GiNaC