A patch for add::coeff() and new clifford_max_label()
Dear All, Working on further patches to clifford.cpp I found that such a simple program #include <iostream> #include <ginac/ginac.h> using namespace std; using namespace GiNaC; #include <stdexcept> int main(){ realsymbol a("a"); varidx nu(symbol("nu", "\\nu"), 2), mu(symbol("mu", "\\mu"), 2); try { ex e = a* dirac_gamma(nu) + 2*dirac_gamma(nu); cout << e.collect(dirac_gamma(nu)) << endl; // -> (2+a)*gamma~nu e = dirac_gamma(mu)* dirac_gamma(nu) + dirac_gamma(nu); cout << e.collect(dirac_gamma(nu)) << endl; // Raised exception: // add::eval(): sum of non-commutative objects has non-zero numeric term } catch (exception &p) { cerr << "Got problem: " << p.what() << endl; } } raises an exception in the second output statement. I think everyone would like it can produce something like "(ONE+gamma~mu)*gamma~nu" instead. And first output should remain "(2+a)*gamma~nu" and not became "(2*ONE+a*ONE)*gamma~nu" on the other hand. I am not sure that is the best way to achieve it but propose a patch which seems to solve it. It introduces one more utility function clifford_max_label() at clifford.cpp, which returns the maximal representation label of non-commutative Clifford object in the expression. Besides that there is also a patch to add::coeff() in add.cpp. I wrote it in mind to reduce unwonted overheads for expression without any Clifford numbers. After this issue will be rectified some more patches to clifford.cpp will follow. Best wishes, Vladimir -- Vladimir V. Kisil email: kisilv@maths.leeds.ac.uk -- www: http://maths.leeds.ac.uk/~kisilv/ Index: ginac/clifford.h =================================================================== RCS file: /home/cvs/GiNaC/ginac/clifford.h,v retrieving revision 1.56 diff -r1.56 clifford.h 285a286,292
/** Returns the maximal representation label of a clifford object * if e contains at least one, otherwise returns -1 * * @param e Expression to be processed * @ignore_ONE defines if clifford_ONE should be ignored in the search*/ ex clifford_max_label(const ex & e, bool ignore_ONE = false);
Index: ginac/clifford.cpp =================================================================== RCS file: /home/cvs/GiNaC/ginac/clifford.cpp,v retrieving revision 1.86 diff -r1.86 clifford.cpp 1035a1084,1100
ex clifford_max_label(const ex & e, bool ignore_ONE) { if (is_a<clifford>(e)) if (ignore_ONE && is_a<diracone>(e.op(0))) return _ex_1; else return ex_to<clifford>(e).get_representation_label(); else { ex rl = _ex_1; for (size_t i=0; i < e.nops(); i++) rl = ex_to<numeric>(rl - clifford_max_label(e.op(i), ignore_ONE)).is_positive()? rl : clifford_max_label(e.op(i), ignore_ONE); return rl; } } Index: ginac/add.cpp =================================================================== RCS file: /home/cvs/GiNaC/ginac/add.cpp,v retrieving revision 1.76 diff -r1.76 add.cpp 31a32 #include "clifford.h" 293a295,297 std::auto_ptr<epvector> coeffseq_cliff(new epvector); ex rl = clifford_max_label(s); bool do_clifford = (!rl.is_equal(_ex_1)), nonscalar = false; 299c303,310 < if (!restcoeff.is_zero())
if (!restcoeff.is_zero()) { if (do_clifford) if (clifford_max_label(restcoeff).is_equal(_ex_1)) coeffseq_cliff->push_back(combine_ex_with_coeff_to_pair(mul(restcoeff,dirac_ONE(ex_to<numeric>(rl).to_int())), i->coeff)); else { coeffseq_cliff->push_back(combine_ex_with_coeff_to_pair(restcoeff, i->coeff)); nonscalar = true; }
300a312
}
304c316,317 < return (new add(coeffseq, n==0 ? overall_coeff : _ex0))->setflag(status_flags::dynallocated); ---
return (new add(nonscalar? coeffseq_cliff : coeffseq, n==0 ? overall_coeff : _ex0))->setflag(status_flags::dynallocated);
Hi Kisil, On Wed, Apr 20, 2005 at 11:07:29AM +0100, Vladimir Kisil wrote:
raises an exception in the second output statement. I think everyone would like it can produce something like "(ONE+gamma~mu)*gamma~nu" instead. And first output should remain "(2+a)*gamma~nu" and not became "(2*ONE+a*ONE)*gamma~nu" on the other hand.
thanks for your patch. It seems to work fine. Nevertheless, I am a bit reluctant to put it into CVS, because there seems to be a bigger issue behind the problem: 'collect' basically uses 'coeff' to do its job. But 'coeff' doesn't care about non-commuting objects. This leads to behavior like the following: ... ex e = dirac_gamma(nu) * dirac_gamma(mu) + dirac_gamma(sigma) * dirac_gamma(mu); cout << e.collect(dirac_gamma(mu)) << endl; gives (gamma~sigma+gamma~nu)*gamma~mu which is fine, but ... ex e = dirac_gamma(mu) * dirac_gamma(nu) + dirac_gamma(mu) * dirac_gamma(sigma); cout << e.collect(dirac_gamma(mu)) << endl; gives (gamma~sigma+gamma~nu)*gamma~mu-gamma~nu*gamma~mu-gamma~sigma*gamma~mu+gamma~mu*gamma~sigma+gamma~mu*gamma~nu which is not wrong but strange. With non-commuting objects one probably needs to distinguish between a left coefficient and a right coefficient. Fixing this problem likely changes the way the ONE-Problem can be addressed. So, at the moment I am not sure what to do. Regards, Jens
Dear Jens, I see the point. To add a puzzle I play with the following test: #include <iostream> #include <ginac/ginac.h> using namespace std; using namespace GiNaC; #include <stdexcept> int main(){ realsymbol a("a"), b("b"), c("c"), d("d"); varidx nu(symbol("nu", "\\nu"), 2), mu(symbol("mu", "\\mu"), 2); try { //cout << latex; ex e = 2*a*b*c*dirac_gamma(nu,1) + 2*dirac_gamma(nu,1)*b*b*a*a; cout << e.collect(dirac_gamma(nu,1)) << endl; // -> gamma~nu*(2*c*b*a+2*b^2*a^2) [OK] e = a*d*dirac_gamma(mu,2)*b*b*a*a* dirac_gamma(nu,1) + b*c*dirac_gamma(nu,1)*d*a+ sin(a)*c*dirac_gamma(nu,1)*abs(b)*b; cout << e.collect(dirac_gamma(nu,1)) << endl; //-> (c*ONE*b*sin(a)*abs(b)+d*c*ONE*b*a+d*gamma~mu*b^2*a^3)*gamma~nu [OK] e = d*dirac_gamma(mu,2) * dirac_gamma(nu,1) + a*b*c*dirac_gamma(nu,1); cout << e.collect(dirac_gamma(nu,1)) << endl; //-> gamma~nu*c*b*a+gamma~nu*(d*gamma~mu+c*ONE*b*a)-gamma~nu*c*ONE*b*a [Ugly!] // Just change representation label of dirac_gamma(mu,1) make iprovement e = d*dirac_gamma(mu,1) * dirac_gamma(nu,1) + a*b*c*dirac_gamma(nu,1); cout << e.collect(dirac_gamma(nu,1)) << endl; //-> (d*gamma~mu+c*ONE*b*a)*gamma~nu [Again OK!??] } catch (exception &p) { cerr << "Got problem: " << p.what() << endl; } } In the third out the expression gamma~nu*c*b*a-gamma~nu*c*ONE*b*a is not simplified to 0 automatically, although this disappear in the fourth output if all Clifford entries has the same representation label. Just in case my previous patch may be of some use I would like to give another version of it. The function clifford_max_label() returns now "char" not "ex" which seems more natural. I think that function may be useful in clifford.cpp anyway. By the way, considering commutativity-noncommutative in GiNaC I realised that all symbols are assumed commutative. Is it difficult to add a new class of noncommutative symbols? Best wishes, Vladimir -- Vladimir V. Kisil email: kisilv@maths.leeds.ac.uk -- www: http://maths.leeds.ac.uk/~kisilv/ Index: ginac/clifford.h =================================================================== RCS file: /home/cvs/GiNaC/ginac/clifford.h,v retrieving revision 1.57 diff -r1.57 clifford.h 287a288,294
/** Returns the maximal representation label of a clifford object * if e contains at least one, otherwise returns -1 * * @param e Expression to be processed * @ignore_ONE defines if clifford_ONE should be ignored in the search*/ char clifford_max_label(const ex & e, bool ignore_ONE = false); Index: ginac/clifford.cpp =================================================================== RCS file: /home/cvs/GiNaC/ginac/clifford.cpp,v retrieving revision 1.88 diff -r1.88 clifford.cpp 1035a1084,1099 char clifford_max_label(const ex & e, bool ignore_ONE) { if (is_a<clifford>(e)) if (ignore_ONE && is_a<diracone>(e.op(0))) return -1; else return ex_to<clifford>(e).get_representation_label(); else { char rl = -1; for (size_t i=0; i < e.nops(); i++) rl = (rl > clifford_max_label(e.op(i), ignore_ONE)) ? rl : clifford_max_label(e.op(i), ignore_ONE); return rl; } }
Index: ginac/add.cpp =================================================================== RCS file: /home/cvs/GiNaC/ginac/add.cpp,v retrieving revision 1.76 diff -r1.76 add.cpp 31a32
#include "clifford.h" 293a295,297 std::auto_ptr<epvector> coeffseq_cliff(new epvector); char rl = clifford_max_label(s); bool do_clifford = (rl != -1), nonscalar = false; 299c303,310 < if (!restcoeff.is_zero())
if (!restcoeff.is_zero()) { if (do_clifford) if (clifford_max_label(restcoeff) == -1) coeffseq_cliff->push_back(combine_ex_with_coeff_to_pair(ncmul(restcoeff, dirac_ONE(rl)), i->coeff)); else { coeffseq_cliff->push_back(combine_ex_with_coeff_to_pair(restcoeff, i->coeff)); nonscalar = true; }
300a312
}
304c316,317 < return (new add(coeffseq, n==0 ? overall_coeff : _ex0))->setflag(status_flags::dynallocated); ---
return (new add(nonscalar? coeffseq_cliff : coeffseq, n==0 ? overall_coeff : _ex0))->setflag(status_flags::dynallocated);
Hi! On Wed, Apr 20, 2005 at 07:11:40PM +0100, Vladimir Kisil wrote:
Is it difficult to add a new class of noncommutative symbols?
You can create noncommutative symbols with this constructor: symbol(const std::string & initname, unsigned rt, unsigned rtt, unsigned domain = domain::complex); "rt" and "rtt" are return_type and return_type_tinfo, respectively. Bye, Christian -- / Physics is an algorithm \/ http://www.uni-mainz.de/~bauec002/
participants (3)
-
Christian Bauer
-
Jens Vollinga
-
Vladimir Kisil