3 * Implementation of GiNaC's color objects.
4 * No real implementation yet, to be done. */
7 * GiNaC Copyright (C) 1999 Johannes Gutenberg University Mainz, Germany
9 * This program is free software; you can redistribute it and/or modify
10 * it under the terms of the GNU General Public License as published by
11 * the Free Software Foundation; either version 2 of the License, or
12 * (at your option) any later version.
14 * This program is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 * GNU General Public License for more details.
19 * You should have received a copy of the GNU General Public License
20 * along with this program; if not, write to the Free Software
21 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
35 #include "relational.h"
38 #ifndef NO_GINAC_NAMESPACE
40 #endif // ndef NO_GINAC_NAMESPACE
43 // default constructor, destructor, copy constructor assignment operator and helpers
48 color::color() : type(invalid), representation_label(0)
50 debugmsg("color default constructor",LOGLEVEL_CONSTRUCT);
51 tinfo_key=TINFO_color;
56 debugmsg("color destructor",LOGLEVEL_DESTRUCT);
60 color::color(color const & other)
62 debugmsg("color copy constructor",LOGLEVEL_CONSTRUCT);
66 color const & color::operator=(color const & other)
68 debugmsg("color operator=",LOGLEVEL_ASSIGNMENT);
78 void color::copy(color const & other)
82 representation_label=other.representation_label;
85 void color::destroy(bool call_parent)
88 indexed::destroy(call_parent);
98 color::color(color_types const t, unsigned const rl) : type(t), representation_label(rl)
100 debugmsg("color constructor from color_types,unsigned",LOGLEVEL_CONSTRUCT);
101 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
102 tinfo_key=TINFO_color;
103 GINAC_ASSERT(all_of_type_coloridx());
106 color::color(color_types const t, ex const & i1, unsigned const rl)
107 : indexed(i1), type(t), representation_label(rl)
109 debugmsg("color constructor from color_types,ex,unsigned",LOGLEVEL_CONSTRUCT);
110 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
111 tinfo_key=TINFO_color;
112 GINAC_ASSERT(all_of_type_coloridx());
115 color::color(color_types const t, ex const & i1, ex const & i2, unsigned const rl)
116 : indexed(i1,i2), type(t), representation_label(rl)
118 debugmsg("color constructor from color_types,ex,ex,unsigned",LOGLEVEL_CONSTRUCT);
119 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
120 tinfo_key=TINFO_color;
121 GINAC_ASSERT(all_of_type_coloridx());
124 color::color(color_types const t, ex const & i1, ex const & i2, ex const & i3,
125 unsigned const rl) : indexed(i1,i2,i3), type(t), representation_label(rl)
127 debugmsg("color constructor from color_types,ex,ex,ex,unsigned",LOGLEVEL_CONSTRUCT);
128 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
129 tinfo_key=TINFO_color;
130 GINAC_ASSERT(all_of_type_coloridx());
133 color::color(color_types const t, exvector const & iv, unsigned const rl)
134 : indexed(iv), type(t), representation_label(rl)
136 debugmsg("color constructor from color_types,exvector,unsigned",LOGLEVEL_CONSTRUCT);
137 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
138 tinfo_key=TINFO_color;
139 GINAC_ASSERT(all_of_type_coloridx());
142 color::color(color_types const t, exvector * ivp, unsigned const rl)
143 : indexed(ivp), type(t), representation_label(rl)
145 debugmsg("color constructor from color_types,exvector *,unsigned",LOGLEVEL_CONSTRUCT);
146 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
147 tinfo_key=TINFO_color;
148 GINAC_ASSERT(all_of_type_coloridx());
152 // functions overriding virtual functions from bases classes
157 basic * color::duplicate() const
159 debugmsg("color duplicate",LOGLEVEL_DUPLICATE);
160 return new color(*this);
163 void color::printraw(ostream & os) const
165 debugmsg("color printraw",LOGLEVEL_PRINT);
166 os << "color(type=" << (unsigned)type
167 << ",representation_label=" << representation_label
170 os << ",hash=" << hashvalue << ",flags=" << flags << ")";
173 void color::printtree(ostream & os, unsigned indent) const
175 debugmsg("color printtree",LOGLEVEL_PRINT);
176 os << string(indent,' ') << "color object: "
177 << "type=" << (unsigned)type
178 << ",representation_label=" << representation_label << ", ";
179 os << seq.size() << " indices" << endl;
180 printtreeindices(os,indent);
181 os << string(indent,' ') << "hash=" << hashvalue
182 << " (0x" << hex << hashvalue << dec << ")"
183 << ", flags=" << flags << endl;
186 void color::print(ostream & os, unsigned upper_precedence) const
188 debugmsg("color print",LOGLEVEL_PRINT);
192 if (representation_label!=0) {
193 os << "^(" << representation_label << ")";
210 os << "INVALID_COLOR_OBJECT";
216 void color::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const
218 debugmsg("color print csrc",LOGLEVEL_PRINT);
219 print(os,upper_precedence);
222 bool color::info(unsigned inf) const
224 return indexed::info(inf);
227 #define CMPINDICES(A,B,C) ((idx1.get_value()==(A))&&(idx2.get_value()==(B))&&(idx3.get_value()==(C)))
229 ex color::eval(int level) const
231 // canonicalize indices
233 bool antisymmetric=false;
237 antisymmetric=true; // no break here!
242 int sig=canonicalize_indices(iv,antisymmetric);
244 // something has changed while sorting indices, more evaluations later
245 if (sig==0) return exZERO();
246 return ex(sig)*color(type,iv,representation_label);
251 // nothing to canonicalize
258 GINAC_ASSERT(seq.size()==2);
259 coloridx const & idx1=ex_to_coloridx(seq[0]);
260 coloridx const & idx2=ex_to_coloridx(seq[1]);
262 // check for delta8_{a,a} where a is a symbolic index, replace by 8
263 if ((idx1.is_symbolic())&&(idx1.is_equal_same_type(idx2))) {
264 return ex(COLOR_EIGHT);
267 // check for delta8_{a,b} where a and b are numeric indices, replace by 0 or 1
268 if ((!idx1.is_symbolic())&&(!idx2.is_symbolic())) {
269 if ((idx1.get_value()!=idx2.get_value())) {
278 // check for d_{a,a,c} (=0) when a is symbolic
280 GINAC_ASSERT(seq.size()==3);
281 coloridx const & idx1=ex_to_coloridx(seq[0]);
282 coloridx const & idx2=ex_to_coloridx(seq[1]);
283 coloridx const & idx3=ex_to_coloridx(seq[2]);
285 if (idx1.is_equal_same_type(idx2) && idx1.is_symbolic()) {
287 } else if (idx2.is_equal_same_type(idx3) && idx2.is_symbolic()) {
291 // check for three numeric indices
292 if (!(idx1.is_symbolic()||idx2.is_symbolic()||idx3.is_symbolic())) {
293 GINAC_ASSERT(idx1.get_value()<=idx2.get_value());
294 GINAC_ASSERT(idx2.get_value()<=idx3.get_value());
295 if (CMPINDICES(1,4,6)||CMPINDICES(1,5,7)||CMPINDICES(2,5,6)||
296 CMPINDICES(3,4,4)||CMPINDICES(3,5,5)) {
298 } else if (CMPINDICES(2,4,7)||CMPINDICES(3,6,6)||CMPINDICES(3,7,7)) {
300 } else if (CMPINDICES(1,1,8)||CMPINDICES(2,2,8)||CMPINDICES(3,3,8)) {
301 return 1/sqrt(numeric(3));
302 } else if (CMPINDICES(8,8,8)) {
303 return -1/sqrt(numeric(3));
304 } else if (CMPINDICES(4,4,8)||CMPINDICES(5,5,8)||CMPINDICES(6,6,8)||CMPINDICES(7,7,8)) {
305 return -1/(2*sqrt(numeric(3)));
313 GINAC_ASSERT(seq.size()==3);
314 coloridx const & idx1=ex_to_coloridx(seq[0]);
315 coloridx const & idx2=ex_to_coloridx(seq[1]);
316 coloridx const & idx3=ex_to_coloridx(seq[2]);
318 // check for three numeric indices
319 if (!(idx1.is_symbolic()||idx2.is_symbolic()||idx3.is_symbolic())) {
320 GINAC_ASSERT(idx1.get_value()<=idx2.get_value());
321 GINAC_ASSERT(idx2.get_value()<=idx3.get_value());
322 if (CMPINDICES(1,2,3)) {
324 } else if (CMPINDICES(1,4,7)||CMPINDICES(2,4,6)||
325 CMPINDICES(2,5,7)||CMPINDICES(3,4,5)) {
327 } else if (CMPINDICES(1,5,6)||CMPINDICES(3,6,7)) {
329 } else if (CMPINDICES(4,5,8)||CMPINDICES(6,7,8)) {
330 return sqrt(numeric(3))/2;
331 } else if (CMPINDICES(8,8,8)) {
332 return -1/sqrt(numeric(3));
333 } else if (CMPINDICES(4,4,8)||CMPINDICES(5,5,8)||CMPINDICES(6,6,8)||CMPINDICES(7,7,8)) {
334 return -1/(2*sqrt(numeric(3)));
341 // nothing to evaluate
350 int color::compare_same_type(basic const & other) const
352 GINAC_ASSERT(other.tinfo() == TINFO_color);
353 const color *o = static_cast<const color *>(&other);
355 if (representation_label==o->representation_label) {
356 return indexed::compare_same_type(other);
358 return representation_label < o->representation_label ? -1 : 1;
360 return type < o->type ? -1 : 1;
363 bool color::is_equal_same_type(basic const & other) const
365 GINAC_ASSERT(other.tinfo() == TINFO_color);
366 const color *o = static_cast<const color *>(&other);
367 if (type!=o->type) return false;
368 if (representation_label!=o->representation_label) return false;
369 return indexed::is_equal_same_type(other);
374 ex color::simplify_ncmul(exvector const & v) const
376 // simplifications: contract delta8_{a,b} where possible
377 // sort delta8,f,d,T(rl=0),T(rl=1),...,ONE(rl=0),ONE(rl=1),...
378 // remove superfluous ONEs
380 // contract indices of delta8_{a,b} if they are different and symbolic
382 exvector v_contracted=v;
383 unsigned replacements;
384 bool something_changed=false;
386 exvector::iterator it=v_contracted.begin();
387 while (it!=v_contracted.end()) {
388 // process only delta8 objects
389 if (is_ex_exactly_of_type(*it,color) && (ex_to_color(*it).type==color_delta8)) {
390 color & d8=ex_to_nonconst_color(*it);
391 GINAC_ASSERT(d8.seq.size()==2);
392 coloridx const & first_idx=ex_to_coloridx(d8.seq[0]);
393 coloridx const & second_idx=ex_to_coloridx(d8.seq[1]);
394 // delta8_{a,a} should have been contracted in color::eval()
395 GINAC_ASSERT((!first_idx.is_equal(second_idx))||(!first_idx.is_symbolic()));
396 ex saved_delta8=*it; // save to restore it later
398 // try to contract first index
400 if (first_idx.is_symbolic()) {
401 replacements = subs_index_in_exvector(v_contracted,first_idx,second_idx);
402 if (replacements==1) {
403 // not contracted except in itself, restore delta8 object
406 // a contracted index should occur exactly twice
407 GINAC_ASSERT(replacements==2);
409 something_changed=true;
413 // try second index only if first was not contracted
414 if ((replacements==1)&&(second_idx.is_symbolic())) {
415 // first index not contracted, *it is guaranteed to be the original delta8 object
416 replacements = subs_index_in_exvector(v_contracted,second_idx,first_idx);
417 if (replacements==1) {
418 // not contracted except in itself, restore delta8 object
421 // a contracted index should occur exactly twice
422 GINAC_ASSERT(replacements==2);
424 something_changed=true;
431 if (something_changed) {
432 // do more simplifications later
433 return nonsimplified_ncmul(v_contracted);
436 // there were no indices to contract
437 // sort delta8,f,d,T(rl=0),T(rl=1),...,ONE(rl=0),ONE(rl=1),...,unknown
438 // (if there is at least one unknown object, all Ts will be unknown to not change the order)
443 exvectorvector Tvecs;
444 Tvecs.resize(MAX_REPRESENTATION_LABELS);
445 exvectorvector ONEvecs;
446 ONEvecs.resize(MAX_REPRESENTATION_LABELS);
449 split_color_string_in_parts(v,delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec);
451 // d_{a,k,l} f_{b,k,l}=0 (includes case a=b)
452 if ((dvec.size()>=1)&&(fvec.size()>=1)) {
453 for (exvector::iterator it1=dvec.begin(); it1!=dvec.end(); ++it1) {
454 for (exvector::iterator it2=fvec.begin(); it2!=fvec.end(); ++it2) {
455 GINAC_ASSERT(is_ex_exactly_of_type(*it1,color));
456 GINAC_ASSERT(is_ex_exactly_of_type(*it2,color));
457 color const & col1=ex_to_color(*it1);
458 color const & col2=ex_to_color(*it2);
459 exvector iv_intersect=idx_intersect(col1.seq,col2.seq);
460 if (iv_intersect.size()>=2) return exZERO();
465 // d_{a,k,l} d_{b,k,l}=5/3 delta8_{a,b} (includes case a=b)
466 if (dvec.size()>=2) {
467 for (exvector::iterator it1=dvec.begin(); it1!=dvec.end()-1; ++it1) {
468 for (exvector::iterator it2=it1+1; it2!=dvec.end(); ++it2) {
469 GINAC_ASSERT(is_ex_exactly_of_type(*it1,color));
470 GINAC_ASSERT(is_ex_exactly_of_type(*it2,color));
471 color const & col1=ex_to_color(*it1);
472 color const & col2=ex_to_color(*it2);
473 exvector iv_intersect=idx_intersect(col1.seq,col2.seq);
474 if (iv_intersect.size()>=2) {
475 if (iv_intersect.size()==3) {
476 *it1=numeric(40)/numeric(3);
479 int sig1, sig2; // unimportant, since symmetric
480 ex idx1=permute_free_index_to_front(col1.seq,iv_intersect,false,&sig1);
481 ex idx2=permute_free_index_to_front(col2.seq,iv_intersect,false,&sig2);
482 *it1=numeric(5)/numeric(3)*color(color_delta8,idx1,idx2);
485 return nonsimplified_ncmul(recombine_color_string(
486 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
492 // f_{a,k,l} f_{b,k,l}=3 delta8_{a,b} (includes case a=b)
493 if (fvec.size()>=2) {
494 for (exvector::iterator it1=fvec.begin(); it1!=fvec.end()-1; ++it1) {
495 for (exvector::iterator it2=it1+1; it2!=fvec.end(); ++it2) {
496 GINAC_ASSERT(is_ex_exactly_of_type(*it1,color));
497 GINAC_ASSERT(is_ex_exactly_of_type(*it2,color));
498 color const & col1=ex_to_color(*it1);
499 color const & col2=ex_to_color(*it2);
500 exvector iv_intersect=idx_intersect(col1.seq,col2.seq);
501 if (iv_intersect.size()>=2) {
502 if (iv_intersect.size()==3) {
507 ex idx1=permute_free_index_to_front(col1.seq,iv_intersect,true,&sig1);
508 ex idx2=permute_free_index_to_front(col2.seq,iv_intersect,true,&sig2);
509 *it1=numeric(sig1*sig2*5)/numeric(3)*color(color_delta8,idx1,idx2);
512 return nonsimplified_ncmul(recombine_color_string(
513 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
519 // d_{a,b,c} T_b T_c = 5/6 T_a
520 // f_{a,b,c} T_b T_c = 3/2 I T_a
521 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
522 if ((Tvecs[rl].size()>=2)&&((dvec.size()>=1)||(fvec.size()>=1))) {
523 for (exvector::iterator it1=Tvecs[rl].begin(); it1!=Tvecs[rl].end()-1; ++it1) {
525 GINAC_ASSERT(is_ex_exactly_of_type(*it1,color)&&ex_to_color(*it1).type==color_T);
526 GINAC_ASSERT(is_ex_exactly_of_type(*(it1+1),color)&&ex_to_color(*(it1+1)).type==color_T);
527 iv.push_back(ex_to_color(*it1).seq[0]);
528 iv.push_back(ex_to_color(*(it1+1)).seq[0]);
530 // d_{a,b,c} T_b T_c = 5/6 T_a
531 for (exvector::iterator it2=dvec.begin(); it2!=dvec.end(); ++it2) {
532 GINAC_ASSERT(is_ex_exactly_of_type(*it2,color)&&ex_to_color(*it2).type==color_d);
533 color const & dref=ex_to_color(*it2);
534 exvector iv_intersect=idx_intersect(dref.seq,iv);
535 if (iv_intersect.size()==2) {
536 int sig; // unimportant, since symmetric
537 ex free_idx=permute_free_index_to_front(dref.seq,iv,false,&sig);
538 *it1=color(color_T,free_idx,rl);
539 *(it1+1)=color(color_ONE,rl);
540 *it2=numeric(5)/numeric(6);
541 return nonsimplified_ncmul(recombine_color_string(
542 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
546 // f_{a,b,c} T_b T_c = 3/2 I T_a
547 for (exvector::iterator it2=fvec.begin(); it2!=fvec.end(); ++it2) {
548 GINAC_ASSERT(is_ex_exactly_of_type(*it2,color)&&ex_to_color(*it2).type==color_f);
549 color const & fref=ex_to_color(*it2);
550 exvector iv_intersect=idx_intersect(fref.seq,iv);
551 if (iv_intersect.size()==2) {
553 ex free_idx=permute_free_index_to_front(fref.seq,iv,true,&sig);
554 *it1=color(color_T,free_idx,rl);
555 *(it1+1)=color(color_ONE,rl);
556 *it2=numeric(sig*3)/numeric(2)*I;
557 return nonsimplified_ncmul(recombine_color_string(
558 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
565 // clear all ONEs when there is at least one corresponding color_T
566 // in this representation, retain one ONE otherwise
567 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
568 if (Tvecs[rl].size()!=0) {
570 } else if (ONEvecs[rl].size()!=0) {
572 ONEvecs[rl].push_back(color(color_ONE,rl));
576 // return a sorted vector
577 return simplified_ncmul(recombine_color_string(delta8vec,fvec,dvec,Tvecs,
578 ONEvecs,unknownvec));
581 ex color::thisexprseq(exvector const & v) const
583 return color(type,v,representation_label);
586 ex color::thisexprseq(exvector * vp) const
588 return color(type,vp,representation_label);
591 bool color::all_of_type_coloridx(void) const
593 // used only inside of ASSERTs
594 for (exvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
595 if (!is_ex_of_type(*cit,coloridx)) return false;
601 // virtual functions which can be overridden by derived classes
607 // non-virtual functions in this class
611 // static member variables
620 const color some_color;
621 type_info const & typeid_color=typeid(some_color);
627 color color_ONE(unsigned const rl)
629 return color(color::color_ONE,rl);
632 color color_T(ex const & a, unsigned const rl)
634 return color(color::color_T,a,rl);
637 color color_f(ex const & a, ex const & b, ex const & c)
639 return color(color::color_f,a,b,c);
642 color color_d(ex const & a, ex const & b, ex const & c)
644 return color(color::color_d,a,b,c);
647 ex color_h(ex const & a, ex const & b, ex const & c)
649 return color(color::color_d,a,b,c)+I*color(color::color_f,a,b,c);
652 color color_delta8(ex const & a, ex const & b)
654 return color(color::color_delta8,a,b);
657 void split_color_string_in_parts(exvector const & v, exvector & delta8vec,
658 exvector & fvec, exvector & dvec,
659 exvectorvector & Tvecs,
660 exvectorvector & ONEvecs,
661 exvector & unknownvec)
663 // if not all elements are of type color, put all Ts in unknownvec to
664 // retain the ordering
666 for (exvector::const_iterator cit=v.begin(); cit!=v.end(); ++cit) {
667 if (!is_ex_exactly_of_type(*cit,color)) {
673 for (exvector::const_iterator cit=v.begin(); cit!=v.end(); ++cit) {
674 if (is_ex_exactly_of_type(*cit,color)) {
675 switch (ex_to_color(*cit).type) {
676 case color::color_delta8:
677 delta8vec.push_back(*cit);
680 fvec.push_back(*cit);
683 dvec.push_back(*cit);
686 GINAC_ASSERT(ex_to_color(*cit).representation_label<MAX_REPRESENTATION_LABELS);
688 Tvecs[ex_to_color(*cit).representation_label].push_back(*cit);
690 unknownvec.push_back(*cit);
693 case color::color_ONE:
694 GINAC_ASSERT(ex_to_color(*cit).representation_label<MAX_REPRESENTATION_LABELS);
695 ONEvecs[ex_to_color(*cit).representation_label].push_back(*cit);
698 throw(std::logic_error("invalid type in split_color_string_in_parts()"));
701 unknownvec.push_back(*cit);
706 exvector recombine_color_string(exvector & delta8vec, exvector & fvec,
707 exvector & dvec, exvectorvector & Tvecs,
708 exvectorvector & ONEvecs, exvector & unknownvec)
710 unsigned sz=delta8vec.size()+fvec.size()+dvec.size()+unknownvec.size();
711 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
712 sz += Tvecs[rl].size();
713 sz += ONEvecs[rl].size();
718 append_exvector_to_exvector(v,delta8vec);
719 append_exvector_to_exvector(v,fvec);
720 append_exvector_to_exvector(v,dvec);
721 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
722 append_exvector_to_exvector(v,Tvecs[rl]);
723 append_exvector_to_exvector(v,ONEvecs[rl]);
725 append_exvector_to_exvector(v,unknownvec);
729 ex color_trace_of_one_representation_label(exvector const & v)
732 return numeric(COLOR_THREE);
733 } else if (v.size()==1) {
734 GINAC_ASSERT(is_ex_exactly_of_type(*(v.begin()),color));
738 ex last_element=v1.back();
739 GINAC_ASSERT(is_ex_exactly_of_type(last_element,color));
740 GINAC_ASSERT(ex_to_color(last_element).type==color::color_T);
742 ex next_to_last_element=v1.back();
743 GINAC_ASSERT(is_ex_exactly_of_type(next_to_last_element,color));
744 GINAC_ASSERT(ex_to_color(next_to_last_element).type==color::color_T);
748 ex const & last_index=ex_to_color(last_element).seq[0];
749 ex const & next_to_last_index=ex_to_color(next_to_last_element).seq[0];
750 ex summation_index=coloridx();
752 v2.push_back(color_T(summation_index)); // don't care about the representation_label
754 // FIXME: check this formula for SU(N) with N!=3
755 return numeric(1)/numeric(2*COLOR_THREE)*color_delta8(next_to_last_index,last_index)
756 % color_trace_of_one_representation_label(v1)
757 +numeric(1)/numeric(2)*color_h(next_to_last_index,last_index,summation_index)
758 % color_trace_of_one_representation_label(v2);
760 ex term1=numeric(1)/numeric(2*COLOR_THREE)*color_delta8(next_to_last_index,last_index)
761 % color_trace_of_one_representation_label(v1);
762 cout << "term 1 of trace of " << v.size() << " ts=" << term1 << endl;
763 ex term2=numeric(1)/numeric(2)*color_h(next_to_last_index,last_index,summation_index)
764 % color_trace_of_one_representation_label(v2);
765 cout << "term 2 of trace of " << v.size() << " ts=" << term2 << endl;
770 ex color_trace(exvector const & v, unsigned const rl)
772 GINAC_ASSERT(rl<MAX_REPRESENTATION_LABELS);
775 v_rest.reserve(v.size()+1); // max size if trace is empty
780 exvectorvector Tvecs;
781 Tvecs.resize(MAX_REPRESENTATION_LABELS);
782 exvectorvector ONEvecs;
783 ONEvecs.resize(MAX_REPRESENTATION_LABELS);
786 split_color_string_in_parts(v,delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec);
788 if (unknownvec.size()!=0) {
789 throw(std::invalid_argument("color_trace(): expression must be expanded"));
792 append_exvector_to_exvector(v_rest,delta8vec);
793 append_exvector_to_exvector(v_rest,fvec);
794 append_exvector_to_exvector(v_rest,dvec);
795 for (unsigned i=0; i<MAX_REPRESENTATION_LABELS; ++i) {
797 append_exvector_to_exvector(v_rest,Tvecs[i]);
798 append_exvector_to_exvector(v_rest,ONEvecs[i]);
800 if (Tvecs[i].size()!=0) {
801 v_rest.push_back(color_trace_of_one_representation_label(Tvecs[i]));
802 } else if (ONEvecs[i].size()!=0) {
803 v_rest.push_back(numeric(COLOR_THREE));
805 throw(std::logic_error("color_trace(): representation_label not in color string"));
810 return nonsimplified_ncmul(v_rest);
813 ex simplify_pure_color_string(ex const & e)
815 GINAC_ASSERT(is_ex_exactly_of_type(e,ncmul));
820 exvectorvector Tvecs;
821 Tvecs.resize(MAX_REPRESENTATION_LABELS);
822 exvectorvector ONEvecs;
823 ONEvecs.resize(MAX_REPRESENTATION_LABELS);
826 split_color_string_in_parts(ex_to_ncmul(e).get_factors(),delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec);
828 // search for T_k S T_k (=1/2 Tr(S) - 1/6 S)
829 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
830 if (Tvecs[rl].size()>=2) {
831 for (unsigned i=0; i<Tvecs[rl].size()-1; ++i) {
832 for (unsigned j=i+1; j<Tvecs[rl].size(); ++j) {
833 ex & t1=Tvecs[rl][i];
834 ex & t2=Tvecs[rl][j];
835 GINAC_ASSERT(is_ex_exactly_of_type(t1,color)&&
836 (ex_to_color(t1).type==color::color_T)&&
837 (ex_to_color(t1).seq.size()==1));
838 GINAC_ASSERT(is_ex_exactly_of_type(t2,color)&&
839 (ex_to_color(t2).type==color::color_T)&&
840 (ex_to_color(t2).seq.size()==1));
841 coloridx const & idx1=ex_to_coloridx(ex_to_color(t1).seq[0]);
842 coloridx const & idx2=ex_to_coloridx(ex_to_color(t2).seq[0]);
844 if (idx1.is_equal(idx2) && idx1.is_symbolic()) {
846 for (unsigned k=i+1; k<j; ++k) {
847 S.push_back(Tvecs[rl][k]);
851 ex term1=numeric(-1)/numeric(6)*nonsimplified_ncmul(recombine_color_string(
852 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
853 for (unsigned k=i+1; k<j; ++k) {
854 S.push_back(exONE());
856 t1=color_trace_of_one_representation_label(S);
857 ex term2=numeric(1)/numeric(2)*nonsimplified_ncmul(recombine_color_string(
858 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
859 return simplify_color(term1+term2);
866 // FIXME: higher contractions
871 ex simplify_color(ex const & e)
873 // all simplification is done on expanded objects
874 ex e_expanded=e.expand();
876 // simplification of sum=sum of simplifications
877 if (is_ex_exactly_of_type(e_expanded,add)) {
879 for (int i=0; i<e_expanded.nops(); ++i) {
880 sum += simplify_color(e_expanded.op(i));
885 // simplification of commutative product=commutative product of simplifications
886 if (is_ex_exactly_of_type(e_expanded,mul)) {
888 for (int i=0; i<e_expanded.nops(); ++i) {
889 prod *= simplify_color(e_expanded.op(i));
894 // simplification of noncommutative product: test if everything is color
895 if (is_ex_exactly_of_type(e_expanded,ncmul)) {
897 for (int i=0; i<e_expanded.nops(); ++i) {
898 if (!is_ex_exactly_of_type(e_expanded.op(i),color)) {
904 return simplify_pure_color_string(e_expanded);
908 // cannot do anything
912 ex brute_force_sum_color_indices(ex const & e)
914 exvector iv_all=e.get_indices();
917 // find double symbolic indices
918 if (iv_all.size()<2) return e;
919 for (exvector::const_iterator cit1=iv_all.begin(); cit1!=iv_all.end()-1; ++cit1) {
920 GINAC_ASSERT(is_ex_of_type(*cit1,coloridx));
921 for (exvector::const_iterator cit2=cit1+1; cit2!=iv_all.end(); ++cit2) {
922 GINAC_ASSERT(is_ex_of_type(*cit2,coloridx));
923 if (ex_to_coloridx(*cit1).is_symbolic() &&
924 ex_to_coloridx(*cit1).is_equal(ex_to_coloridx(*cit2))) {
925 iv_double.push_back(*cit1);
932 counter.resize(iv_double.size());
934 for (l=0; unsigned(l)<iv_double.size(); ++l) {
942 for (l=0; unsigned(l)<iv_double.size(); ++l) {
943 term=term.subs(iv_double[l]==coloridx((unsigned)(counter[l])));
944 //iv_double[l].print(cout);
945 //cout << " " << counter[l] << " ";
950 // increment counter[]
951 l=iv_double.size()-1;
952 while ((l>=0)&&((++counter[l])>(int)COLOR_EIGHT)) {
956 if (l<2) { cout << counter[0] << counter[1] << endl; }
963 void append_exvector_to_exvector(exvector & dest, exvector const & source)
965 for (exvector::const_iterator cit=source.begin(); cit!=source.end(); ++cit) {
966 dest.push_back(*cit);
970 #ifndef NO_GINAC_NAMESPACE
972 #endif // ndef NO_GINAC_NAMESPACE