From: Chris Dams Date: Fri, 20 Oct 2006 16:29:33 +0000 (+0000) Subject: Bugs w.r.t. raising/lowering dummy indices when doing simplify_indexed. X-Git-Tag: release_1-4-0~55 X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?p=ginac.git;a=commitdiff_plain;h=f01c43a7b5cffea04d34ba3d30b41781255ad5d4;hp=d4df6322c4790ea932280602fdb584486f6101b6;ds=sidebyside Bugs w.r.t. raising/lowering dummy indices when doing simplify_indexed. --- diff --git a/check/exam_clifford.cpp b/check/exam_clifford.cpp index d6544eda..e8ee3f70 100644 --- a/check/exam_clifford.cpp +++ b/check/exam_clifford.cpp @@ -64,9 +64,14 @@ static unsigned check_equal_simplify_term(const ex & e1, const ex & e2, idx & mu ex e = expand_dummy_sum(normal(simplify_indexed(e1) - e2), true); for (int j=0; j<4; j++) { - ex esub = e.subs(lst(is_a(mu) ? - mu == idx(j, mu.get_dim()), ex_to(mu).toggle_variance() == idx(j, mu.get_dim()) - : mu == idx(j, mu.get_dim()))); + ex esub = e.subs( + is_a(mu) + ? lst ( + mu == idx(j, mu.get_dim()), + ex_to(mu).toggle_variance() == idx(j, mu.get_dim()) + ) + : lst(mu == idx(j, mu.get_dim())) + ); if (!(canonicalize_clifford(esub).is_zero())) { clog << "simplify_indexed(" << e1 << ") - (" << e2 << ") erroneously returned " << canonicalize_clifford(esub) << " instead of 0 for mu=" << j << endl; @@ -312,151 +317,139 @@ static unsigned clifford_check5() * the cases when used indexes have or have not variance. * To this end we recycle the code through the following macros */ -#define CHECK6(IDX,TOGGLE) {IDX v(symbol("v"), 4), nu(symbol("nu"), 4), mu(symbol("mu"), 4), \ - psi(symbol("psi"),4), lam(symbol("lambda"), 4),\ - xi(symbol("xi"), 4), rho(symbol("rho"),4);\ -\ -/* checks general identities and contractions for clifford_unit*/\ - e = dirac_ONE(2) * clifford_unit(mu, A, 2) * dirac_ONE(2);\ - result += check_equal(e, clifford_unit(mu, A, 2));\ -\ - e = clifford_unit(idx(2, 4), A) * clifford_unit(idx(1, 4), A)\ - * clifford_unit(idx(1, 4), A) * clifford_unit(idx(2, 4), A);\ - result += check_equal(e, A(1, 1) * A(2, 2) * dirac_ONE());\ -\ - e = clifford_unit(IDX(2, 4), A) * clifford_unit(IDX(1, 4), A)\ - * clifford_unit(IDX(1, 4), A) * clifford_unit(IDX(2, 4), A);\ - result += check_equal(e, A(1, 1) * A(2, 2) * dirac_ONE());\ -\ - e = clifford_unit(nu, A) * clifford_unit(nu TOGGLE, A);\ - result += check_equal_simplify(e, A.trace() * dirac_ONE());\ -\ - e = clifford_unit(nu, A) * clifford_unit(nu, A);\ - result += check_equal_simplify(e, indexed(A_symm, sy_symm(), nu, nu) * dirac_ONE());\ -\ - e = clifford_unit(nu, A) * clifford_unit(nu TOGGLE, A) * clifford_unit(mu, A);\ - result += check_equal_simplify(e, A.trace() * clifford_unit(mu, A));\ -\ - e = clifford_unit(nu, A) * clifford_unit(mu, A) * clifford_unit(nu TOGGLE, A);\ - \ - result += check_equal_simplify_term(e, 2 * indexed(A_symm, sy_symm(), nu TOGGLE, mu) *clifford_unit(nu, A)-A.trace()*clifford_unit(mu, A), mu);\ -\ - e = clifford_unit(nu, A) * clifford_unit(nu TOGGLE, A)\ - * clifford_unit(mu, A) * clifford_unit(mu TOGGLE, A);\ - result += check_equal_simplify(e, pow(A.trace(), 2) * dirac_ONE());\ -\ - e = clifford_unit(mu, A) * clifford_unit(nu, A)\ - * clifford_unit(nu TOGGLE, A) * clifford_unit(mu TOGGLE, A);\ - result += check_equal_simplify(e, pow(A.trace(), 2) * dirac_ONE());\ -\ - e = clifford_unit(mu, A) * clifford_unit(nu, A)\ - * clifford_unit(mu TOGGLE, A) * clifford_unit(nu TOGGLE, A);\ -\ - result += check_equal_simplify_term2(e, 2*indexed(A_symm, sy_symm(), nu TOGGLE, mu TOGGLE) * clifford_unit(mu, A) * clifford_unit(nu, A) - pow(A.trace(), 2)*dirac_ONE());\ -\ - e = clifford_unit(mu TOGGLE, A) * clifford_unit(nu, A)\ - * clifford_unit(mu, A) * clifford_unit(nu TOGGLE, A);\ -\ - result += check_equal_simplify_term2(e, 2*indexed(A_symm, nu, mu) * clifford_unit(mu TOGGLE, A) * clifford_unit(nu TOGGLE, A) - pow(A.trace(), 2)*dirac_ONE());\ -\ - e = clifford_unit(nu TOGGLE, A) * clifford_unit(rho TOGGLE, A)\ - * clifford_unit(mu, A) * clifford_unit(rho, A) * clifford_unit(nu, A);\ - e = e.simplify_indexed().collect(clifford_unit(mu, A));\ - \ - result += check_equal_simplify_term(e, 4* indexed(A_symm, sy_symm(), nu TOGGLE, rho)*indexed(A_symm, sy_symm(), rho TOGGLE, mu) *clifford_unit(nu, A) \ - - 2*A.trace() * (clifford_unit(rho, A) * indexed(A_symm, sy_symm(), rho TOGGLE, mu) \ - + clifford_unit(nu, A) * indexed(A_symm, sy_symm(), nu TOGGLE, mu)) + pow(A.trace(),2)* clifford_unit(mu, A), mu);\ -\ - e = clifford_unit(nu TOGGLE, A) * clifford_unit(rho, A)\ - * clifford_unit(mu, A) * clifford_unit(rho TOGGLE, A) * clifford_unit(nu, A);\ - e = e.simplify_indexed().collect(clifford_unit(mu, A));\ - \ - result += check_equal_simplify_term(e, 4* indexed(A_symm, sy_symm(), nu TOGGLE, rho)*indexed(A_symm, sy_symm(), rho TOGGLE, mu) *clifford_unit(nu, A) \ - - 2*A.trace() * (clifford_unit(rho, A) * indexed(A_symm, sy_symm(), rho TOGGLE, mu) \ - + clifford_unit(nu, A) * indexed(A_symm, sy_symm(), nu TOGGLE, mu)) + pow(A.trace(),2)* clifford_unit(mu, A), mu);\ -\ - e = clifford_unit(mu, A) * clifford_unit(nu, A) + clifford_unit(nu, A) * clifford_unit(mu, A);\ - result += check_equal(canonicalize_clifford(e), 2*dirac_ONE()*indexed(A_symm, sy_symm(), mu, nu));\ -\ - e = (clifford_unit(mu, A) * clifford_unit(nu, A) * clifford_unit(lam, A)\ - + clifford_unit(nu, A) * clifford_unit(lam, A) * clifford_unit(mu, A)\ - + clifford_unit(lam, A) * clifford_unit(mu, A) * clifford_unit(nu, A)\ - - clifford_unit(nu, A) * clifford_unit(mu, A) * clifford_unit(lam, A)\ - - clifford_unit(lam, A) * clifford_unit(nu, A) * clifford_unit(mu, A)\ - - clifford_unit(mu, A) * clifford_unit(lam, A) * clifford_unit(nu, A)) / 6\ - + indexed(A_symm, sy_symm(), mu, nu) * clifford_unit(lam, A)\ - - indexed(A_symm, sy_symm(), mu, lam) * clifford_unit(nu, A)\ - + indexed(A_symm, sy_symm(), nu, lam) * clifford_unit(mu, A)\ - - clifford_unit(mu, A) * clifford_unit(nu, A) * clifford_unit(lam, A);\ - result += check_equal(canonicalize_clifford(e), 0);\ -\ -/* lst_to_clifford() and clifford_inverse() check*/\ - realsymbol x("x"), y("y"), t("t"), z("z");\ - \ - ex c = clifford_unit(nu, A, 1);\ - e = lst_to_clifford(lst(t, x, y, z), mu, A, 1) * lst_to_clifford(lst(1, 2, 3, 4), c);\ - e1 = clifford_inverse(e);\ - result += check_equal_simplify_term2((e*e1).simplify_indexed(), dirac_ONE(1));\ -\ -/* Moebius map (both forms) checks for symmetric metrics only */\ - matrix M1(2, 2), M2(2, 2);\ - c = clifford_unit(nu, A);\ - \ - e = clifford_moebius_map(0, dirac_ONE(), \ - dirac_ONE(), 0, lst(t, x, y, z), A); \ -/* this is just the inversion*/\ - M1 = 0, dirac_ONE(),\ - dirac_ONE(), 0;\ - e1 = clifford_moebius_map(M1, lst(t, x, y, z), A); \ -/* the inversion again*/\ - result += check_equal_lst(e, e1);\ - \ - e1 = clifford_to_lst(clifford_inverse(lst_to_clifford(lst(t, x, y, z), mu, A)), c);\ - result += check_equal_lst(e, e1);\ - \ - e = clifford_moebius_map(dirac_ONE(), lst_to_clifford(lst(1, 2, 3, 4), nu, A), \ - 0, dirac_ONE(), lst(t, x, y, z), A); \ -/*this is just a shift*/\ - M2 = dirac_ONE(), lst_to_clifford(lst(1, 2, 3, 4), c),\ - 0, dirac_ONE();\ - e1 = clifford_moebius_map(M2, lst(t, x, y, z), c); \ -/* the same shift*/\ - result += check_equal_lst(e, e1);\ - \ - result += check_equal(e, lst(t+1, x+2, y+3, z+4));\ - \ -/* Check the group law for Moebius maps */\ - e = clifford_moebius_map(M1, ex_to(e1), c); \ -/*composition of M1 and M2*/\ - e1 = clifford_moebius_map(M1.mul(M2), lst(t, x, y, z), c); \ -/* the product M1*M2*/\ - result += check_equal_lst(e, e1);} - -static unsigned clifford_check6(const matrix & A) +template unsigned clifford_check6(const matrix &A) { + unsigned result = 0; + matrix A_symm(4,4), A2(4, 4); A_symm = A.add(A.transpose()).mul(half); A2 = A_symm.mul(A_symm); - + + IDX v(symbol("v"), 4), nu(symbol("nu"), 4), mu(symbol("mu"), 4), + psi(symbol("psi"),4), lam(symbol("lambda"), 4), + xi(symbol("xi"), 4), rho(symbol("rho"),4); + ex mu_TOGGLE = is_a(mu) ? ex_to(mu).toggle_variance() : mu; + ex nu_TOGGLE = is_a(nu) ? ex_to(nu).toggle_variance() : nu; + ex rho_TOGGLE + = is_a(rho) ? ex_to(rho).toggle_variance() : rho; + ex e, e1; - int result = 0; + +/* checks general identities and contractions for clifford_unit*/ + e = dirac_ONE(2) * clifford_unit(mu, A, 2) * dirac_ONE(2); + result += check_equal(e, clifford_unit(mu, A, 2)); + + e = clifford_unit(idx(2, 4), A) * clifford_unit(idx(1, 4), A) + * clifford_unit(idx(1, 4), A) * clifford_unit(idx(2, 4), A); + result += check_equal(e, A(1, 1) * A(2, 2) * dirac_ONE()); + + e = clifford_unit(IDX(2, 4), A) * clifford_unit(IDX(1, 4), A) + * clifford_unit(IDX(1, 4), A) * clifford_unit(IDX(2, 4), A); + result += check_equal(e, A(1, 1) * A(2, 2) * dirac_ONE()); + + e = clifford_unit(nu, A) * clifford_unit(nu_TOGGLE, A); + result += check_equal_simplify(e, A.trace() * dirac_ONE()); + + e = clifford_unit(nu, A) * clifford_unit(nu, A); + result += check_equal_simplify(e, indexed(A_symm, sy_symm(), nu, nu) * dirac_ONE()); + + e = clifford_unit(nu, A) * clifford_unit(nu_TOGGLE, A) * clifford_unit(mu, A); + result += check_equal_simplify(e, A.trace() * clifford_unit(mu, A)); + + e = clifford_unit(nu, A) * clifford_unit(mu, A) * clifford_unit(nu_TOGGLE, A); - CHECK6(varidx,.toggle_variance()) + result += check_equal_simplify_term(e, 2 * indexed(A_symm, sy_symm(), nu_TOGGLE, mu) *clifford_unit(nu, A)-A.trace()*clifford_unit(mu, A), mu); - return result; -} + e = clifford_unit(nu, A) * clifford_unit(nu_TOGGLE, A) + * clifford_unit(mu, A) * clifford_unit(mu_TOGGLE, A); + result += check_equal_simplify(e, pow(A.trace(), 2) * dirac_ONE()); -static unsigned clifford_check6a(const matrix & A) -{ - matrix A_symm(4,4), A2(4, 4); - A_symm = A.add(A.transpose()).mul(half); - A2 = A_symm.mul(A_symm); + e = clifford_unit(mu, A) * clifford_unit(nu, A) + * clifford_unit(nu_TOGGLE, A) * clifford_unit(mu_TOGGLE, A); + result += check_equal_simplify(e, pow(A.trace(), 2) * dirac_ONE()); + + e = clifford_unit(mu, A) * clifford_unit(nu, A) + * clifford_unit(mu_TOGGLE, A) * clifford_unit(nu_TOGGLE, A); + + result += check_equal_simplify_term2(e, 2*indexed(A_symm, sy_symm(), nu_TOGGLE, mu_TOGGLE) * clifford_unit(mu, A) * clifford_unit(nu, A) - pow(A.trace(), 2)*dirac_ONE()); + + e = clifford_unit(mu_TOGGLE, A) * clifford_unit(nu, A) + * clifford_unit(mu, A) * clifford_unit(nu_TOGGLE, A); + + result += check_equal_simplify_term2(e, 2*indexed(A_symm, nu, mu) * clifford_unit(mu_TOGGLE, A) * clifford_unit(nu_TOGGLE, A) - pow(A.trace(), 2)*dirac_ONE()); + + e = clifford_unit(nu_TOGGLE, A) * clifford_unit(rho_TOGGLE, A) + * clifford_unit(mu, A) * clifford_unit(rho, A) * clifford_unit(nu, A); + e = e.simplify_indexed().collect(clifford_unit(mu, A)); - ex e, e1; - int result = 0; + result += check_equal_simplify_term(e, 4* indexed(A_symm, sy_symm(), nu_TOGGLE, rho)*indexed(A_symm, sy_symm(), rho_TOGGLE, mu) *clifford_unit(nu, A) + - 2*A.trace() * (clifford_unit(rho, A) * indexed(A_symm, sy_symm(), rho_TOGGLE, mu) + + clifford_unit(nu, A) * indexed(A_symm, sy_symm(), nu_TOGGLE, mu)) + pow(A.trace(),2)* clifford_unit(mu, A), mu); - CHECK6(idx,) + e = clifford_unit(nu_TOGGLE, A) * clifford_unit(rho, A) + * clifford_unit(mu, A) * clifford_unit(rho_TOGGLE, A) * clifford_unit(nu, A); + e = e.simplify_indexed().collect(clifford_unit(mu, A)); + + result += check_equal_simplify_term(e, 4* indexed(A_symm, sy_symm(), nu_TOGGLE, rho)*indexed(A_symm, sy_symm(), rho_TOGGLE, mu) *clifford_unit(nu, A) + - 2*A.trace() * (clifford_unit(rho, A) * indexed(A_symm, sy_symm(), rho_TOGGLE, mu) + + clifford_unit(nu, A) * indexed(A_symm, sy_symm(), nu_TOGGLE, mu)) + pow(A.trace(),2)* clifford_unit(mu, A), mu); + + e = clifford_unit(mu, A) * clifford_unit(nu, A) + clifford_unit(nu, A) * clifford_unit(mu, A); + result += check_equal(canonicalize_clifford(e), 2*dirac_ONE()*indexed(A_symm, sy_symm(), mu, nu)); + + e = (clifford_unit(mu, A) * clifford_unit(nu, A) * clifford_unit(lam, A) + + clifford_unit(nu, A) * clifford_unit(lam, A) * clifford_unit(mu, A) + + clifford_unit(lam, A) * clifford_unit(mu, A) * clifford_unit(nu, A) + - clifford_unit(nu, A) * clifford_unit(mu, A) * clifford_unit(lam, A) + - clifford_unit(lam, A) * clifford_unit(nu, A) * clifford_unit(mu, A) + - clifford_unit(mu, A) * clifford_unit(lam, A) * clifford_unit(nu, A)) / 6 + + indexed(A_symm, sy_symm(), mu, nu) * clifford_unit(lam, A) + - indexed(A_symm, sy_symm(), mu, lam) * clifford_unit(nu, A) + + indexed(A_symm, sy_symm(), nu, lam) * clifford_unit(mu, A) + - clifford_unit(mu, A) * clifford_unit(nu, A) * clifford_unit(lam, A); + result += check_equal(canonicalize_clifford(e), 0); +/* lst_to_clifford() and clifford_inverse() check*/ + realsymbol x("x"), y("y"), t("t"), z("z"); + + ex c = clifford_unit(nu, A, 1); + e = lst_to_clifford(lst(t, x, y, z), mu, A, 1) * lst_to_clifford(lst(1, 2, 3, 4), c); + e1 = clifford_inverse(e); + result += check_equal_simplify_term2((e*e1).simplify_indexed(), dirac_ONE(1)); + +/* Moebius map (both forms) checks for symmetric metrics only */ + matrix M1(2, 2), M2(2, 2); + c = clifford_unit(nu, A); + + e = clifford_moebius_map(0, dirac_ONE(), + dirac_ONE(), 0, lst(t, x, y, z), A); +/* this is just the inversion*/ + M1 = 0, dirac_ONE(), + dirac_ONE(), 0; + e1 = clifford_moebius_map(M1, lst(t, x, y, z), A); +/* the inversion again*/ + result += check_equal_lst(e, e1); + + e1 = clifford_to_lst(clifford_inverse(lst_to_clifford(lst(t, x, y, z), mu, A)), c); + result += check_equal_lst(e, e1); + + e = clifford_moebius_map(dirac_ONE(), lst_to_clifford(lst(1, 2, 3, 4), nu, A), + 0, dirac_ONE(), lst(t, x, y, z), A); +/*this is just a shift*/ + M2 = dirac_ONE(), lst_to_clifford(lst(1, 2, 3, 4), c), + 0, dirac_ONE(); + e1 = clifford_moebius_map(M2, lst(t, x, y, z), c); +/* the same shift*/ + result += check_equal_lst(e, e1); + + result += check_equal(e, lst(t+1, x+2, y+3, z+4)); + +/* Check the group law for Moebius maps */ + e = clifford_moebius_map(M1, ex_to(e1), c); +/*composition of M1 and M2*/ + e1 = clifford_moebius_map(M1.mul(M2), lst(t, x, y, z), c); +/* the product M1*M2*/ + result += check_equal_lst(e, e1); return result; } @@ -545,45 +538,46 @@ unsigned exam_clifford() result += clifford_check5(); cout << '.' << flush; // anticommuting, symmetric examples - result += clifford_check6(ex_to(diag_matrix(lst(-1, 1, 1, 1))))+clifford_check6a(ex_to(diag_matrix(lst(-1, 1, 1, 1))));; cout << '.' << flush; - result += clifford_check6(ex_to(diag_matrix(lst(-1, -1, -1, -1))))+clifford_check6a(ex_to(diag_matrix(lst(-1, -1, -1, -1))));; cout << '.' << flush; - result += clifford_check6(ex_to(diag_matrix(lst(-1, 1, 1, -1))))+clifford_check6a(ex_to(diag_matrix(lst(-1, 1, 1, -1))));; cout << '.' << flush; - result += clifford_check6(ex_to(diag_matrix(lst(-1, 0, 1, -1))))+clifford_check6a(ex_to(diag_matrix(lst(-1, 0, 1, -1))));; cout << '.' << flush; - result += clifford_check6(ex_to(diag_matrix(lst(-3, 0, 2, -1))))+clifford_check6a(ex_to(diag_matrix(lst(-3, 0, 2, -1))));; cout << '.' << flush; + result += clifford_check6(ex_to(diag_matrix(lst(-1, 1, 1, 1)))); + result += clifford_check6(ex_to(diag_matrix(lst(-1, 1, 1, 1))));; cout << '.' << flush; + result += clifford_check6(ex_to(diag_matrix(lst(-1, -1, -1, -1))))+clifford_check6(ex_to(diag_matrix(lst(-1, -1, -1, -1))));; cout << '.' << flush; + result += clifford_check6(ex_to(diag_matrix(lst(-1, 1, 1, -1))))+clifford_check6(ex_to(diag_matrix(lst(-1, 1, 1, -1))));; cout << '.' << flush; + result += clifford_check6(ex_to(diag_matrix(lst(-1, 0, 1, -1))))+clifford_check6(ex_to(diag_matrix(lst(-1, 0, 1, -1))));; cout << '.' << flush; + result += clifford_check6(ex_to(diag_matrix(lst(-3, 0, 2, -1))))+clifford_check6(ex_to(diag_matrix(lst(-3, 0, 2, -1))));; cout << '.' << flush; realsymbol s("s"), t("t"); // symbolic entries in matric - result += clifford_check6(ex_to(diag_matrix(lst(-1, 1, s, t))))+clifford_check6a(ex_to(diag_matrix(lst(-1, 1, s, t))));; cout << '.' << flush; + result += clifford_check6(ex_to(diag_matrix(lst(-1, 1, s, t))))+clifford_check6(ex_to(diag_matrix(lst(-1, 1, s, t))));; cout << '.' << flush; matrix A(4, 4); A = 1, 0, 0, 0, // anticommuting, not symmetric, Tr=0 0, -1, 0, 0, 0, 0, 0, -1, 0, 0, 1, 0; - result += clifford_check6(A)+clifford_check6a(A);; cout << '.' << flush; + result += clifford_check6(A)+clifford_check6(A);; cout << '.' << flush; A = 1, 0, 0, 0, // anticommuting, not symmetric, Tr=2 0, 1, 0, 0, 0, 0, 0, -1, 0, 0, 1, 0; - result += clifford_check6(A)+clifford_check6a(A);; cout << '.' << flush; + result += clifford_check6(A)+clifford_check6(A);; cout << '.' << flush; A = 1, 0, 0, 0, // not anticommuting, symmetric, Tr=0 0, -1, 0, 0, 0, 0, 0, -1, 0, 0, -1, 0; - result += clifford_check6(A)+clifford_check6a(A);; cout << '.' << flush; + result += clifford_check6(A)+clifford_check6(A);; cout << '.' << flush; A = 1, 0, 0, 0, // not anticommuting, symmetric, Tr=2 0, 1, 0, 0, 0, 0, 0, -1, 0, 0, -1, 0; - result += clifford_check6(A)+clifford_check6a(A);; cout << '.' << flush; + result += clifford_check6(A)+clifford_check6(A);; cout << '.' << flush; A = 1, 1, 0, 0, // not anticommuting, not symmetric, Tr=4 0, 1, 1, 0, 0, 0, 1, 1, 0, 0, 0, 1; - result += clifford_check6(A)+clifford_check6a(A);; cout << '.' << flush; + result += clifford_check6(A)+clifford_check6(A);; cout << '.' << flush; symbol dim("D"); result += clifford_check7(minkmetric(), dim); cout << '.' << flush; diff --git a/check/exam_paranoia.cpp b/check/exam_paranoia.cpp index 7c82e486..0ab64910 100644 --- a/check/exam_paranoia.cpp +++ b/check/exam_paranoia.cpp @@ -437,6 +437,33 @@ static unsigned exam_paranoia16() return result; } +// Bug in reposition_dummy_indices() could result in correct expression +// turned into one with inconsistent indices. Fixed on Aug 29, 2006 +static unsigned exam_paranoia17() +{ + varidx mu1(symbol("mu1"), 4); + varidx mu2(symbol("mu2"), 4); + varidx mu3(symbol("mu3"), 4); + varidx mu4(symbol("mu4"), 4); + varidx mu5(symbol("mu5"), 4); + varidx mu6(symbol("mu6"), 4); + + exvector ev2; + ev2.push_back(mu3.toggle_variance()); + ev2.push_back(mu6); + ev2.push_back(mu5.toggle_variance()); + ev2.push_back(mu6.toggle_variance()); + ev2.push_back(mu5); + ev2.push_back(mu3); + // notice: all indices are contracted ... + + ex test_cycl = indexed(symbol("A"), sy_cycl(), ev2); + test_cycl = test_cycl.simplify_indexed(); + // ... so there should be zero free indices in the end. + return test_cycl.get_free_indices().size(); +} + + unsigned exam_paranoia() { unsigned result = 0; @@ -460,6 +487,7 @@ unsigned exam_paranoia() result += exam_paranoia14(); cout << '.' << flush; result += exam_paranoia15(); cout << '.' << flush; result += exam_paranoia16(); cout << '.' << flush; + result += exam_paranoia17(); cout << '.' << flush; if (!result) { cout << " passed " << endl; diff --git a/ginac/indexed.cpp b/ginac/indexed.cpp index d772b3f4..445c273e 100644 --- a/ginac/indexed.cpp +++ b/ginac/indexed.cpp @@ -639,60 +639,60 @@ bool reposition_dummy_indices(ex & e, exvector & variant_dummy_indices, exvector { bool something_changed = false; - // Canonicalize wrt indices that are dummies within e. I.e., their - // symbol occurs twice in an index of e. This is only done if there - // is a cyclic symmetry because in that case it may happen that after - // raising/lowering an index the indices get reshuffled by ::eval in - // such a way that the iterators no longer point to the right objects. - if (ex_to(ex_to(e).get_symmetry()).has_cyclic()) { - // Get dummy pairs of varidxes within the indexed object in e. - exvector local_var_dummies; - local_var_dummies.reserve(e.nops()/2); - for (size_t i=1; i(e.op(i))) - continue; - for (size_t j=i+1; jop(0)) { - variant_dummy_indices.erase(k); - break; - } + // Find dummy symbols that occur twice in the same indexed object. + exvector local_var_dummies; + local_var_dummies.reserve(e.nops()/2); + for (size_t i=1; i(e.op(i))) + continue; + for (size_t j=i+1; jop(0)) { + variant_dummy_indices.erase(k); + break; } - break; } + break; } } - // Try all posibilities of raising/lowering and keep the least one in - // the sense of ex_is_less. - ex optimal_e = e; - size_t numpossibs = 1 << local_var_dummies.size(); - for (size_t i=0; i(curr_idx).toggle_variance(); - m[curr_idx] = curr_toggle; - m[curr_toggle] = curr_idx; - } - try_e = e.subs(m, subs_options::no_pattern); - } - if(ex_is_less()(try_e, optimal_e)) - { optimal_e = try_e; - something_changed = true; + } + + // In the case where a dummy symbol occurs twice in the same indexed object + // we try all posibilities of raising/lowering and keep the least one in + // the sense of ex_is_less. + ex optimal_e = e; + size_t numpossibs = 1 << local_var_dummies.size(); + for (size_t i=0; i(curr_idx).toggle_variance(); + m[curr_idx] = curr_toggle; + m[curr_toggle] = curr_idx; } + try_e = e.subs(m, subs_options::no_pattern); + } + if(ex_is_less()(try_e, optimal_e)) + { optimal_e = try_e; + something_changed = true; } - e = optimal_e; } + e = optimal_e; + + if (!is_a(e)) + return true; + + exvector seq = ex_to(e).seq; // If a dummy index is encountered for the first time in the // product, pull it up, otherwise, pull it down - exvector::const_iterator it2, it2start, it2end; - for (it2start = ex_to(e).seq.begin(), it2end = ex_to(e).seq.end(), it2 = it2start + 1; it2 != it2end; ++it2) { + for (exvector::iterator it2 = seq.begin()+1, it2end = seq.end(); + it2 != it2end; ++it2) { if (!is_exactly_a(*it2)) continue; @@ -700,14 +700,20 @@ bool reposition_dummy_indices(ex & e, exvector & variant_dummy_indices, exvector for (vit = variant_dummy_indices.begin(), vitend = variant_dummy_indices.end(); vit != vitend; ++vit) { if (it2->op(0).is_equal(vit->op(0))) { if (ex_to(*it2).is_covariant()) { - e = e.subs(lst( - *it2 == ex_to(*it2).toggle_variance(), - ex_to(*it2).toggle_variance() == *it2 - ), subs_options::no_pattern); + /* + * N.B. we don't want to use + * + * e = e.subs(lst( + * *it2 == ex_to(*it2).toggle_variance(), + * ex_to(*it2).toggle_variance() == *it2 + * ), subs_options::no_pattern); + * + * since this can trigger non-trivial repositioning of indices, + * e.g. due to non-trivial symmetry properties of e, thus + * invalidating iterators + */ + *it2 = ex_to(*it2).toggle_variance(); something_changed = true; - it2 = ex_to(e).seq.begin() + (it2 - it2start); - it2start = ex_to(e).seq.begin(); - it2end = ex_to(e).seq.end(); } moved_indices.push_back(*vit); variant_dummy_indices.erase(vit); @@ -718,11 +724,8 @@ bool reposition_dummy_indices(ex & e, exvector & variant_dummy_indices, exvector for (vit = moved_indices.begin(), vitend = moved_indices.end(); vit != vitend; ++vit) { if (it2->op(0).is_equal(vit->op(0))) { if (ex_to(*it2).is_contravariant()) { - e = e.subs(*it2 == ex_to(*it2).toggle_variance(), subs_options::no_pattern); + *it2 = ex_to(*it2).toggle_variance(); something_changed = true; - it2 = ex_to(e).seq.begin() + (it2 - it2start); - it2start = ex_to(e).seq.begin(); - it2end = ex_to(e).seq.end(); } goto next_index; } @@ -731,6 +734,9 @@ bool reposition_dummy_indices(ex & e, exvector & variant_dummy_indices, exvector next_index: ; } + if (something_changed) + e = ex_to(e).thiscontainer(seq); + return something_changed; }