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";
}