[GiNaC-list] Substitution of indexed expressions

Z S baltcode at gmail.com
Mon Nov 16 00:37:29 CET 2009

I am new to ginac, and I want to simplify some expected values of indexed
expressions. As a toy example, I wanted to simplify
E[\sum_{i,j} A.i.j ^ 2]
where it is known that E[A.i.j * A.i.j] = 1.
The answer should be N + E[ \sum{i,j} A.i.j \sum_{k != i,l != j} A.k.l].

I tried to use the following program but it doesnt work. Any help would be
greatly appreciated as I need to implement similar but more involved
expressions. Thanks!

#include <iostream>
#include <ginac/ginac.h>
using namespace std;
using namespace GiNaC;

int main()
{
symbol N("N");
symbol i_sym("i"), j_sym("j");
idx i(i_sym, N), j(j_sym, N);
symbol A("A");
ex e = indexed(A,i,j) * indexed(A,i,j);
cout << "e = "<< e<< endl;
cout << "Free indices = " << exprseq(e.get_free_indices()) << "\n";

// I want to substitute A.i.j * A.i.j = 1, for any i,j
// The following gives a run time error
//ex e3 = subs(e,indexed(A,wild(),wild()) * indexed(A,wild(),wild())
== 1);

// The following doesnt substitute/simplify anything

idx i0(wild(),N), i1(wild(),N);
ex e3 = subs(e,indexed(A,i0,i1) * indexed(A,i0,i1) == 1);
cout << "e3 = " << e3 << endl;
cout << "Free indices = " << exprseq(e3.get_free_indices()) << "\n";
}
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.cebix.net/pipermail/ginac-list/attachments/20091115/5dea2627/attachment.html>