* Implementation of GiNaC's indexed expressions. */
/*
- * GiNaC Copyright (C) 1999-2001 Johannes Gutenberg University Mainz, Germany
+ * GiNaC Copyright (C) 1999-2002 Johannes Gutenberg University Mainz, Germany
*
* This program is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
*/
+#include <iostream>
#include <stdexcept>
-#include <algorithm>
#include "indexed.h"
#include "idx.h"
#include "print.h"
#include "archive.h"
#include "utils.h"
-#include "debugmsg.h"
namespace GiNaC {
GINAC_IMPLEMENT_REGISTERED_CLASS(indexed, exprseq)
//////////
-// default constructor, destructor, copy constructor assignment operator and helpers
+// default ctor, dtor, copy ctor, assignment operator and helpers
//////////
indexed::indexed() : symtree(sy_none())
{
- debugmsg("indexed default constructor", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
}
indexed::indexed(const ex & b) : inherited(b), symtree(sy_none())
{
- debugmsg("indexed constructor from ex", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
validate();
}
indexed::indexed(const ex & b, const ex & i1) : inherited(b, i1), symtree(sy_none())
{
- debugmsg("indexed constructor from ex,ex", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
validate();
}
indexed::indexed(const ex & b, const ex & i1, const ex & i2) : inherited(b, i1, i2), symtree(sy_none())
{
- debugmsg("indexed constructor from ex,ex,ex", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
validate();
}
indexed::indexed(const ex & b, const ex & i1, const ex & i2, const ex & i3) : inherited(b, i1, i2, i3), symtree(sy_none())
{
- debugmsg("indexed constructor from ex,ex,ex,ex", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
validate();
}
indexed::indexed(const ex & b, const ex & i1, const ex & i2, const ex & i3, const ex & i4) : inherited(b, i1, i2, i3, i4), symtree(sy_none())
{
- debugmsg("indexed constructor from ex,ex,ex,ex,ex", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
validate();
}
indexed::indexed(const ex & b, const symmetry & symm, const ex & i1, const ex & i2) : inherited(b, i1, i2), symtree(symm)
{
- debugmsg("indexed constructor from ex,symmetry,ex,ex", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
validate();
}
indexed::indexed(const ex & b, const symmetry & symm, const ex & i1, const ex & i2, const ex & i3) : inherited(b, i1, i2, i3), symtree(symm)
{
- debugmsg("indexed constructor from ex,symmetry,ex,ex,ex", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
validate();
}
indexed::indexed(const ex & b, const symmetry & symm, const ex & i1, const ex & i2, const ex & i3, const ex & i4) : inherited(b, i1, i2, i3, i4), symtree(symm)
{
- debugmsg("indexed constructor from ex,symmetry,ex,ex,ex,ex", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
validate();
}
indexed::indexed(const ex & b, const exvector & v) : inherited(b), symtree(sy_none())
{
- debugmsg("indexed constructor from ex,exvector", LOGLEVEL_CONSTRUCT);
seq.insert(seq.end(), v.begin(), v.end());
tinfo_key = TINFO_indexed;
validate();
indexed::indexed(const ex & b, const symmetry & symm, const exvector & v) : inherited(b), symtree(symm)
{
- debugmsg("indexed constructor from ex,symmetry,exvector", LOGLEVEL_CONSTRUCT);
seq.insert(seq.end(), v.begin(), v.end());
tinfo_key = TINFO_indexed;
validate();
indexed::indexed(const symmetry & symm, const exprseq & es) : inherited(es), symtree(symm)
{
- debugmsg("indexed constructor from symmetry,exprseq", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
}
indexed::indexed(const symmetry & symm, const exvector & v, bool discardable) : inherited(v, discardable), symtree(symm)
{
- debugmsg("indexed constructor from symmetry,exvector", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
}
indexed::indexed(const symmetry & symm, exvector * vp) : inherited(vp), symtree(symm)
{
- debugmsg("indexed constructor from symmetry,exvector *", LOGLEVEL_CONSTRUCT);
tinfo_key = TINFO_indexed;
}
indexed::indexed(const archive_node &n, const lst &sym_lst) : inherited(n, sym_lst)
{
- debugmsg("indexed constructor from archive_node", LOGLEVEL_CONSTRUCT);
if (!n.find_ex("symmetry", symtree, sym_lst)) {
// GiNaC versions <= 0.9.0 had an unsigned "symmetry" property
unsigned symm = 0;
symtree = sy_none();
break;
}
- ex_to_nonconst_symmetry(symtree).validate(seq.size() - 1);
+ const_cast<symmetry &>(ex_to<symmetry>(symtree)).validate(seq.size() - 1);
}
}
void indexed::print(const print_context & c, unsigned level) const
{
- debugmsg("indexed print", LOGLEVEL_PRINT);
GINAC_ASSERT(seq.size() > 0);
if (is_of_type(c, print_tree)) {
int indexed::compare_same_type(const basic & other) const
{
- GINAC_ASSERT(is_of_type(other, indexed));
+ GINAC_ASSERT(is_a<indexed>(other));
return inherited::compare_same_type(other);
}
// If the base object is 0, the whole object is 0
if (base.is_zero())
- return _ex0();
+ return _ex0;
// 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)) {
// Canonicalize indices according to the symmetry properties
if (seq.size() > 2) {
exvector v = seq;
- GINAC_ASSERT(is_ex_exactly_of_type(symtree, symmetry));
+ GINAC_ASSERT(is_exactly_a<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)
- return _ex0();
+ return _ex0;
return ex(sig) * thisexprseq(v);
}
}
// Let the class of the base object perform additional evaluations
- return base.bp->eval_indexed(*this);
-}
-
-int indexed::degree(const ex & s) const
-{
- return is_equal(*s.bp) ? 1 : 0;
-}
-
-int indexed::ldegree(const ex & s) const
-{
- return is_equal(*s.bp) ? 1 : 0;
-}
-
-ex indexed::coeff(const ex & s, int n) const
-{
- if (is_equal(*s.bp))
- return n==1 ? _ex1() : _ex0();
- else
- return n==0 ? ex(*this) : _ex0();
+ return ex_to<basic>(base).eval_indexed(*this);
}
ex indexed::thisexprseq(const exvector & v) const
// expand_indexed expands (a+b).i -> a.i + b.i
const ex & base = seq[0];
- ex sum = _ex0();
+ ex sum = _ex0;
for (unsigned i=0; i<base.nops(); i++) {
exvector s = seq;
s[0] = base.op(i);
while (it != itend) {
bool cur_covariant = (is_ex_of_type(*it, varidx) ? ex_to<varidx>(*it).is_covariant() : true);
- if (first || cur_covariant != covariant) {
+ if (first || cur_covariant != covariant) { // Variance changed
+ // The empty {} prevents indices from ending up on top of each other
if (!first)
- c.s << "}";
+ c.s << "}{}";
covariant = cur_covariant;
if (covariant)
c.s << "_{";
if (!symtree.is_zero()) {
if (!is_ex_exactly_of_type(symtree, symmetry))
throw(std::invalid_argument("symmetry of indexed object must be of type symmetry"));
- ex_to_nonconst_symmetry(symtree).validate(seq.size() - 1);
+ const_cast<symmetry &>(ex_to<symmetry>(symtree)).validate(seq.size() - 1);
}
}
* @see ex::diff */
ex indexed::derivative(const symbol & s) const
{
- return _ex0();
+ return _ex0;
}
//////////
}
it++;
}
- shaker_sort(global_dummy_indices.begin(), global_dummy_indices.end(), ex_is_less(), ex_swap());
// If this is the first set of local indices, do nothing
if (old_global_size == 0)
shaker_sort(local_syms.begin(), local_syms.end(), ex_is_less(), ex_swap());
for (unsigned i=0; i<global_size; i++)
global_syms.push_back(global_dummy_indices[i].op(0));
+ shaker_sort(global_syms.begin(), global_syms.end(), ex_is_less(), ex_swap());
// Remove common indices
exlist local_uniq, global_uniq;
if (is_ex_exactly_of_type(e, power)) {
// We only get called for simple squares, split a^2 -> a*a
- GINAC_ASSERT(e.op(1).is_equal(_ex2()));
+ GINAC_ASSERT(e.op(1).is_equal(_ex2));
v.push_back(e.op(0));
v.push_back(e.op(0));
} else {
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())) {
+ if (is_ex_exactly_of_type(f, power) && f.op(1).is_equal(_ex2)) {
v.push_back(f.op(0));
v.push_back(f.op(0));
} else if (is_ex_exactly_of_type(f, ncmul)) {
if (free.empty()) {
if (sp.is_defined(*it1, *it2)) {
*it1 = sp.evaluate(*it1, *it2);
- *it2 = _ex1();
+ *it2 = _ex1;
goto contraction_done;
}
}
// Try to contract the first one with the second one
- contracted = it1->op(0).bp->contract_with(it1, it2, v);
+ contracted = ex_to<basic>(it1->op(0)).contract_with(it1, it2, v);
if (!contracted) {
// That didn't work; maybe the second object knows how to
// contract itself with the first one
- contracted = it2->op(0).bp->contract_with(it2, it1, v);
+ contracted = ex_to<basic>(it2->op(0)).contract_with(it2, it1, v);
}
if (contracted) {
contraction_done:
dummy_syms.append(local_dummy_indices[i].op(0));
if (r.symmetrize(dummy_syms).is_zero()) {
free_indices.clear();
- return _ex0();
+ return _ex0;
}
}
// 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 ex_to<basic>(r.op(0).op(0)).scalar_mul_indexed(r.op(0), ex_to<numeric>(r.op(1)));
else
return r;
}
// free indices in each term
if (is_ex_exactly_of_type(e_expanded, add)) {
bool first = true;
- ex sum = _ex0();
+ ex sum = _ex0;
free_indices.clear();
for (unsigned i=0; i<e_expanded.nops(); i++) {
if (!indices_consistent(free_indices, free_indices_of_term))
throw (std::runtime_error("simplify_indexed: inconsistent indices in sum"));
if (is_ex_of_type(sum, indexed) && is_ex_of_type(term, indexed))
- sum = sum.op(0).bp->add_indexed(sum, term);
+ sum = ex_to<basic>(sum.op(0)).add_indexed(sum, term);
else
sum += term;
}
// Simplification of products
if (is_ex_exactly_of_type(e_expanded, mul)
|| is_ex_exactly_of_type(e_expanded, ncmul)
- || (is_ex_exactly_of_type(e_expanded, power) && is_ex_of_type(e_expanded.op(0), indexed) && e_expanded.op(1).is_equal(_ex2())))
+ || (is_ex_exactly_of_type(e_expanded, power) && is_ex_of_type(e_expanded.op(0), indexed) && e_expanded.op(1).is_equal(_ex2)))
return simplify_indexed_product(e_expanded, free_indices, dummy_indices, sp);
// Cannot do anything