]> www.ginac.de Git - ginac.git/blobdiff - ginac/indexed.cpp
- insert class fderivative.
[ginac.git] / ginac / indexed.cpp
index fd96eeb9cab71bec23146e65551b6ce26a341e59..c6a84e31ac70e6346091651bd0eb61f1a70fef5d 100644 (file)
@@ -187,7 +187,7 @@ void indexed::archive(archive_node &n) const
 DEFAULT_UNARCHIVE(indexed)
 
 //////////
-// functions overriding virtual functions from bases classes
+// functions overriding virtual functions from base classes
 //////////
 
 void indexed::print(const print_context & c, unsigned level) const
@@ -201,7 +201,6 @@ void indexed::print(const print_context & c, unsigned level) const
                    << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
                    << ", " << seq.size()-1 << " indices"
                    << ", symmetry=" << symtree << std::endl;
-               c.s << std::endl;
                unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
                seq[0].print(c, level + delta_indent);
                printindices(c, level + delta_indent);
@@ -235,7 +234,7 @@ bool indexed::info(unsigned inf) const
 
 struct idx_is_not : public std::binary_function<ex, unsigned, bool> {
        bool operator() (const ex & e, unsigned inf) const {
-               return !(ex_to_idx(e).get_value().info(inf));
+               return !(ex_to<idx>(e).get_value().info(inf));
        }
 };
 
@@ -259,7 +258,7 @@ ex indexed::eval(int level) const
 {
        // First evaluate children, then we will end up here again
        if (level > 1)
-               return indexed(ex_to_symmetry(symtree), evalchildren(level));
+               return indexed(ex_to<symmetry>(symtree), evalchildren(level));
 
        const ex &base = seq[0];
 
@@ -270,7 +269,7 @@ ex indexed::eval(int level) const
        // If the base object is a product, pull out the numeric factor
        if (is_ex_exactly_of_type(base, mul) && is_ex_exactly_of_type(base.op(base.nops() - 1), numeric)) {
                exvector v(seq);
-               ex f = ex_to_numeric(base.op(base.nops() - 1));
+               ex f = ex_to<numeric>(base.op(base.nops() - 1));
                v[0] = seq[0] / f;
                return f * thisexprseq(v);
        }
@@ -279,7 +278,7 @@ ex indexed::eval(int level) const
        if (seq.size() > 2) {
                exvector v = seq;
                GINAC_ASSERT(is_ex_exactly_of_type(symtree, symmetry));
-               int sig = canonicalize(v.begin() + 1, ex_to_symmetry(symtree));
+               int sig = canonicalize(v.begin() + 1, ex_to<symmetry>(symtree));
                if (sig != INT_MAX) {
                        // Something has changed while sorting indices, more evaluations later
                        if (sig == 0)
@@ -312,12 +311,12 @@ ex indexed::coeff(const ex & s, int n) const
 
 ex indexed::thisexprseq(const exvector & v) const
 {
-       return indexed(ex_to_symmetry(symtree), v);
+       return indexed(ex_to<symmetry>(symtree), v);
 }
 
 ex indexed::thisexprseq(exvector * vp) const
 {
-       return indexed(ex_to_symmetry(symtree), vp);
+       return indexed(ex_to<symmetry>(symtree), vp);
 }
 
 ex indexed::expand(unsigned options) const
@@ -363,7 +362,7 @@ void indexed::printindices(const print_context & c, unsigned level) const
                        bool covariant = true;
 
                        while (it != itend) {
-                               bool cur_covariant = (is_ex_of_type(*it, varidx) ? ex_to_varidx(*it).is_covariant() : true);
+                               bool cur_covariant = (is_ex_of_type(*it, varidx) ? ex_to<varidx>(*it).is_covariant() : true);
                                if (first || cur_covariant != covariant) {
                                        if (!first)
                                                c.s << "}";
@@ -411,6 +410,14 @@ void indexed::validate(void) const
        }
 }
 
+/** Implementation of ex::diff() for an indexed object always returns 0.
+ *
+ *  @see ex::diff */
+ex indexed::derivative(const symbol & s) const
+{
+       return _ex0();
+}
+
 //////////
 // global functions
 //////////
@@ -527,8 +534,8 @@ exvector power::get_free_indices(void) const
  *    by the function */
 static ex rename_dummy_indices(const ex & e, exvector & global_dummy_indices, exvector & local_dummy_indices)
 {
-       int global_size = global_dummy_indices.size(),
-           local_size = local_dummy_indices.size();
+       unsigned global_size = global_dummy_indices.size(),
+                local_size = local_dummy_indices.size();
 
        // Any local dummy indices at all?
        if (local_size == 0)
@@ -557,7 +564,8 @@ static ex rename_dummy_indices(const ex & e, exvector & global_dummy_indices, ex
        for (unsigned i=0; i<local_size; i++) {
                ex loc_sym = local_dummy_indices[i].op(0);
                ex glob_sym = global_dummy_indices[i].op(0);
-               if (!loc_sym.is_equal(glob_sym)) {
+               if (!loc_sym.is_equal(glob_sym)
+                && ex_to<idx>(local_dummy_indices[i]).get_dim().is_equal(ex_to<idx>(global_dummy_indices[i]).get_dim())) {
                        all_equal = false;
                        local_syms.append(loc_sym);
                        global_syms.append(glob_sym);
@@ -587,7 +595,7 @@ ex simplify_indexed_product(const ex & e, exvector & free_indices, exvector & du
                v.push_back(e.op(0));
                v.push_back(e.op(0));
        } else {
-               for (int i=0; i<e.nops(); i++) {
+               for (unsigned i=0; i<e.nops(); i++) {
                        ex f = e.op(i);
                        if (is_ex_exactly_of_type(f, power) && f.op(1).is_equal(_ex2())) {
                                v.push_back(f.op(0));
@@ -595,7 +603,7 @@ ex simplify_indexed_product(const ex & e, exvector & free_indices, exvector & du
                        } else if (is_ex_exactly_of_type(f, ncmul)) {
                                // Noncommutative factor found, split it as well
                                non_commutative = true; // everything becomes noncommutative, ncmul will sort out the commutative factors later
-                               for (int j=0; j<f.nops(); j++)
+                               for (unsigned j=0; j<f.nops(); j++)
                                        v.push_back(f.op(j));
                        } else
                                v.push_back(f);
@@ -617,7 +625,7 @@ try_again:
                // Indexed factor found, get free indices and look for contraction
                // candidates
                exvector free1, dummy1;
-               find_free_and_dummy(ex_to_indexed(*it1).seq.begin() + 1, ex_to_indexed(*it1).seq.end(), free1, dummy1);
+               find_free_and_dummy(ex_to<indexed>(*it1).seq.begin() + 1, ex_to<indexed>(*it1).seq.end(), free1, dummy1);
 
                exvector::iterator it2;
                for (it2 = it1 + 1; it2 != itend; it2++) {
@@ -630,18 +638,19 @@ try_again:
                        // Find free indices of second factor and merge them with free
                        // indices of first factor
                        exvector un;
-                       find_free_and_dummy(ex_to_indexed(*it2).seq.begin() + 1, ex_to_indexed(*it2).seq.end(), un, dummy1);
+                       find_free_and_dummy(ex_to<indexed>(*it2).seq.begin() + 1, ex_to<indexed>(*it2).seq.end(), un, dummy1);
                        un.insert(un.end(), free1.begin(), free1.end());
 
                        // Check whether the two factors share dummy indices
                        exvector free, dummy;
                        find_free_and_dummy(un, free, dummy);
-                       if (dummy.size() == 0)
+                       unsigned num_dummies = dummy.size();
+                       if (num_dummies == 0)
                                continue;
 
                        // At least one dummy index, is it a defined scalar product?
                        bool contracted = false;
-                       if (free.size() == 0) {
+                       if (free.empty()) {
                                if (sp.is_defined(*it1, *it2)) {
                                        *it1 = sp.evaluate(*it1, *it2);
                                        *it2 = _ex1();
@@ -650,13 +659,13 @@ try_again:
                        }
 
                        // Contraction of symmetric with antisymmetric object is zero
-                       if (dummy.size() > 1
-                        && ex_to_symmetry(ex_to_indexed(*it1).symtree).has_symmetry()
-                        && ex_to_symmetry(ex_to_indexed(*it2).symtree).has_symmetry()) {
+                       if (num_dummies > 1
+                        && ex_to<symmetry>(ex_to<indexed>(*it1).symtree).has_symmetry()
+                        && ex_to<symmetry>(ex_to<indexed>(*it2).symtree).has_symmetry()) {
 
                                // Check all pairs of dummy indices
-                               for (unsigned idx1=0; idx1<dummy.size()-1; idx1++) {
-                                       for (unsigned idx2=idx1+1; idx2<dummy.size(); idx2++) {
+                               for (unsigned idx1=0; idx1<num_dummies-1; idx1++) {
+                                       for (unsigned idx2=idx1+1; idx2<num_dummies; idx2++) {
 
                                                // Try and swap the index pair and check whether the
                                                // relative sign changed
@@ -713,7 +722,7 @@ contraction_done:
                exvector free_indices_of_factor;
                if (is_ex_of_type(*it1, indexed)) {
                        exvector dummy_indices_of_factor;
-                       find_free_and_dummy(ex_to_indexed(*it1).seq.begin() + 1, ex_to_indexed(*it1).seq.end(), free_indices_of_factor, dummy_indices_of_factor);
+                       find_free_and_dummy(ex_to<indexed>(*it1).seq.begin() + 1, ex_to<indexed>(*it1).seq.end(), free_indices_of_factor, dummy_indices_of_factor);
                        individual_dummy_indices.insert(individual_dummy_indices.end(), dummy_indices_of_factor.begin(), dummy_indices_of_factor.end());
                } else
                        free_indices_of_factor = it1->get_free_indices();
@@ -736,7 +745,7 @@ contraction_done:
        // Product of indexed object with a scalar?
        if (is_ex_exactly_of_type(r, mul) && r.nops() == 2
         && is_ex_exactly_of_type(r.op(1), numeric) && is_ex_of_type(r.op(0), indexed))
-               return r.op(0).op(0).bp->scalar_mul_indexed(r.op(0), ex_to_numeric(r.op(1)));
+               return r.op(0).op(0).bp->scalar_mul_indexed(r.op(0), ex_to<numeric>(r.op(1)));
        else
                return r;
 }
@@ -750,7 +759,7 @@ ex simplify_indexed(const ex & e, exvector & free_indices, exvector & dummy_indi
        // Simplification of single indexed object: just find the free indices
        // and perform dummy index renaming
        if (is_ex_of_type(e_expanded, indexed)) {
-               const indexed &i = ex_to_indexed(e_expanded);
+               const indexed &i = ex_to<indexed>(e_expanded);
                exvector local_dummy_indices;
                find_free_and_dummy(i.seq.begin() + 1, i.seq.end(), free_indices, local_dummy_indices);
                return rename_dummy_indices(e_expanded, dummy_indices, local_dummy_indices);
@@ -881,10 +890,12 @@ ex scalar_products::evaluate(const ex & v1, const ex & v2) const
 void scalar_products::debugprint(void) const
 {
        std::cerr << "map size=" << spm.size() << std::endl;
-       for (spmap::const_iterator cit=spm.begin(); cit!=spm.end(); ++cit) {
-               const spmapkey & k = cit->first;
+       spmap::const_iterator i = spm.begin(), end = spm.end();
+       while (i != end) {
+               const spmapkey & k = i->first;
                std::cerr << "item key=(" << k.first << "," << k.second;
-               std::cerr << "), value=" << cit->second << std::endl;
+               std::cerr << "), value=" << i->second << std::endl;
+               ++i;
        }
 }