3 * Implementation of GiNaC's color objects.
4 * No real implementation yet, to be done. */
15 // default constructor, destructor, copy constructor assignment operator and helpers
20 color::color() : type(invalid), representation_label(0)
22 debugmsg("color default constructor",LOGLEVEL_CONSTRUCT);
23 tinfo_key=TINFO_COLOR;
28 debugmsg("color destructor",LOGLEVEL_DESTRUCT);
32 color::color(color const & other)
34 debugmsg("color copy constructor",LOGLEVEL_CONSTRUCT);
38 color const & color::operator=(color const & other)
40 debugmsg("color operator=",LOGLEVEL_ASSIGNMENT);
50 void color::copy(color const & other)
54 representation_label=other.representation_label;
57 void color::destroy(bool call_parent)
60 indexed::destroy(call_parent);
70 color::color(color_types const t, unsigned const rl) : type(t), representation_label(rl)
72 debugmsg("color constructor from color_types,unsigned",LOGLEVEL_CONSTRUCT);
73 ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
74 tinfo_key=TINFO_COLOR;
75 ASSERT(all_of_type_coloridx());
78 color::color(color_types const t, ex const & i1, unsigned const rl)
79 : indexed(i1), type(t), representation_label(rl)
81 debugmsg("color constructor from color_types,ex,unsigned",LOGLEVEL_CONSTRUCT);
82 ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
83 tinfo_key=TINFO_COLOR;
84 ASSERT(all_of_type_coloridx());
87 color::color(color_types const t, ex const & i1, ex const & i2, unsigned const rl)
88 : indexed(i1,i2), type(t), representation_label(rl)
90 debugmsg("color constructor from color_types,ex,ex,unsigned",LOGLEVEL_CONSTRUCT);
91 ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
92 tinfo_key=TINFO_COLOR;
93 ASSERT(all_of_type_coloridx());
96 color::color(color_types const t, ex const & i1, ex const & i2, ex const & i3,
97 unsigned const rl) : indexed(i1,i2,i3), type(t), representation_label(rl)
99 debugmsg("color constructor from color_types,ex,ex,ex,unsigned",LOGLEVEL_CONSTRUCT);
100 ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
101 tinfo_key=TINFO_COLOR;
102 ASSERT(all_of_type_coloridx());
105 color::color(color_types const t, exvector const & iv, unsigned const rl)
106 : indexed(iv), type(t), representation_label(rl)
108 debugmsg("color constructor from color_types,exvector,unsigned",LOGLEVEL_CONSTRUCT);
109 ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
110 tinfo_key=TINFO_COLOR;
111 ASSERT(all_of_type_coloridx());
114 color::color(color_types const t, exvector * ivp, unsigned const rl)
115 : indexed(ivp), type(t), representation_label(rl)
117 debugmsg("color constructor from color_types,exvector *,unsigned",LOGLEVEL_CONSTRUCT);
118 ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
119 tinfo_key=TINFO_COLOR;
120 ASSERT(all_of_type_coloridx());
124 // functions overriding virtual functions from bases classes
129 basic * color::duplicate() const
131 debugmsg("color duplicate",LOGLEVEL_DUPLICATE);
132 return new color(*this);
135 void color::printraw(ostream & os) const
137 debugmsg("color printraw",LOGLEVEL_PRINT);
138 os << "color(type=" << (unsigned)type
139 << ",representation_label=" << representation_label
142 os << ",hash=" << hashvalue << ",flags=" << flags << ")";
145 void color::printtree(ostream & os, unsigned indent) const
147 debugmsg("color printtree",LOGLEVEL_PRINT);
148 os << string(indent,' ') << "color object: "
149 << "type=" << (unsigned)type
150 << ",representation_label=" << representation_label << ", ";
151 os << seq.size() << " indices" << endl;
152 printtreeindices(os,indent);
153 os << string(indent,' ') << "hash=" << hashvalue
154 << " (0x" << hex << hashvalue << dec << ")"
155 << ", flags=" << flags << endl;
158 void color::print(ostream & os, unsigned upper_precedence) const
160 debugmsg("color print",LOGLEVEL_PRINT);
164 if (representation_label!=0) {
165 os << "^(" << representation_label << ")";
182 os << "INVALID_COLOR_OBJECT";
188 void color::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const
190 debugmsg("color print csrc",LOGLEVEL_PRINT);
191 print(os,upper_precedence);
194 bool color::info(unsigned inf) const
196 return indexed::info(inf);
199 #define CMPINDICES(A,B,C) ((idx1.get_value()==(A))&&(idx2.get_value()==(B))&&(idx3.get_value()==(C)))
201 ex color::eval(int level) const
203 // canonicalize indices
205 bool antisymmetric=false;
209 antisymmetric=true; // no break here!
214 int sig=canonicalize_indices(iv,antisymmetric);
216 // something has changed while sorting indices, more evaluations later
217 if (sig==0) return exZERO();
218 return ex(sig)*color(type,iv,representation_label);
223 // nothing to canonicalize
230 ASSERT(seq.size()==2);
231 coloridx const & idx1=ex_to_coloridx(seq[0]);
232 coloridx const & idx2=ex_to_coloridx(seq[1]);
234 // check for delta8_{a,a} where a is a symbolic index, replace by 8
235 if ((idx1.is_symbolic())&&(idx1.is_equal_same_type(idx2))) {
236 return ex(COLOR_EIGHT);
239 // check for delta8_{a,b} where a and b are numeric indices, replace by 0 or 1
240 if ((!idx1.is_symbolic())&&(!idx2.is_symbolic())) {
241 if ((idx1.get_value()!=idx2.get_value())) {
250 // check for d_{a,a,c} (=0) when a is symbolic
252 ASSERT(seq.size()==3);
253 coloridx const & idx1=ex_to_coloridx(seq[0]);
254 coloridx const & idx2=ex_to_coloridx(seq[1]);
255 coloridx const & idx3=ex_to_coloridx(seq[2]);
257 if (idx1.is_equal_same_type(idx2) && idx1.is_symbolic()) {
259 } else if (idx2.is_equal_same_type(idx3) && idx2.is_symbolic()) {
263 // check for three numeric indices
264 if (!(idx1.is_symbolic()||idx2.is_symbolic()||idx3.is_symbolic())) {
265 ASSERT(idx1.get_value()<=idx2.get_value());
266 ASSERT(idx2.get_value()<=idx3.get_value());
267 if (CMPINDICES(1,4,6)||CMPINDICES(1,5,7)||CMPINDICES(2,5,6)||
268 CMPINDICES(3,4,4)||CMPINDICES(3,5,5)) {
270 } else if (CMPINDICES(2,4,7)||CMPINDICES(3,6,6)||CMPINDICES(3,7,7)) {
272 } else if (CMPINDICES(1,1,8)||CMPINDICES(2,2,8)||CMPINDICES(3,3,8)) {
273 return 1/sqrt(numeric(3));
274 } else if (CMPINDICES(8,8,8)) {
275 return -1/sqrt(numeric(3));
276 } else if (CMPINDICES(4,4,8)||CMPINDICES(5,5,8)||CMPINDICES(6,6,8)||CMPINDICES(7,7,8)) {
277 return -1/(2*sqrt(numeric(3)));
285 ASSERT(seq.size()==3);
286 coloridx const & idx1=ex_to_coloridx(seq[0]);
287 coloridx const & idx2=ex_to_coloridx(seq[1]);
288 coloridx const & idx3=ex_to_coloridx(seq[2]);
290 // check for three numeric indices
291 if (!(idx1.is_symbolic()||idx2.is_symbolic()||idx3.is_symbolic())) {
292 ASSERT(idx1.get_value()<=idx2.get_value());
293 ASSERT(idx2.get_value()<=idx3.get_value());
294 if (CMPINDICES(1,2,3)) {
296 } else if (CMPINDICES(1,4,7)||CMPINDICES(2,4,6)||
297 CMPINDICES(2,5,7)||CMPINDICES(3,4,5)) {
299 } else if (CMPINDICES(1,5,6)||CMPINDICES(3,6,7)) {
301 } else if (CMPINDICES(4,5,8)||CMPINDICES(6,7,8)) {
302 return sqrt(numeric(3))/2;
303 } else if (CMPINDICES(8,8,8)) {
304 return -1/sqrt(numeric(3));
305 } else if (CMPINDICES(4,4,8)||CMPINDICES(5,5,8)||CMPINDICES(6,6,8)||CMPINDICES(7,7,8)) {
306 return -1/(2*sqrt(numeric(3)));
313 // nothing to evaluate
322 int color::compare_same_type(basic const & other) const
324 ASSERT(other.tinfo() == TINFO_COLOR);
325 const color *o = static_cast<const color *>(&other);
327 if (representation_label==o->representation_label) {
328 return indexed::compare_same_type(other);
330 return representation_label < o->representation_label ? -1 : 1;
332 return type < o->type ? -1 : 1;
335 bool color::is_equal_same_type(basic const & other) const
337 ASSERT(other.tinfo() == TINFO_COLOR);
338 const color *o = static_cast<const color *>(&other);
339 if (type!=o->type) return false;
340 if (representation_label!=o->representation_label) return false;
341 return indexed::is_equal_same_type(other);
346 ex color::simplify_ncmul(exvector const & v) const
348 // simplifications: contract delta8_{a,b} where possible
349 // sort delta8,f,d,T(rl=0),T(rl=1),...,ONE(rl=0),ONE(rl=1),...
350 // remove superfluous ONEs
352 // contract indices of delta8_{a,b} if they are different and symbolic
354 exvector v_contracted=v;
355 unsigned replacements;
356 bool something_changed=false;
358 exvector::iterator it=v_contracted.begin();
359 while (it!=v_contracted.end()) {
360 // process only delta8 objects
361 if (is_ex_exactly_of_type(*it,color) && (ex_to_color(*it).type==color_delta8)) {
362 color & d8=ex_to_nonconst_color(*it);
363 ASSERT(d8.seq.size()==2);
364 coloridx const & first_idx=ex_to_coloridx(d8.seq[0]);
365 coloridx const & second_idx=ex_to_coloridx(d8.seq[1]);
366 // delta8_{a,a} should have been contracted in color::eval()
367 ASSERT((!first_idx.is_equal(second_idx))||(!first_idx.is_symbolic()));
368 ex saved_delta8=*it; // save to restore it later
370 // try to contract first index
372 if (first_idx.is_symbolic()) {
373 replacements = subs_index_in_exvector(v_contracted,first_idx,second_idx);
374 if (replacements==1) {
375 // not contracted except in itself, restore delta8 object
378 // a contracted index should occur exactly twice
379 ASSERT(replacements==2);
381 something_changed=true;
385 // try second index only if first was not contracted
386 if ((replacements==1)&&(second_idx.is_symbolic())) {
387 // first index not contracted, *it is guaranteed to be the original delta8 object
388 replacements = subs_index_in_exvector(v_contracted,second_idx,first_idx);
389 if (replacements==1) {
390 // not contracted except in itself, restore delta8 object
393 // a contracted index should occur exactly twice
394 ASSERT(replacements==2);
396 something_changed=true;
403 if (something_changed) {
404 // do more simplifications later
405 return nonsimplified_ncmul(v_contracted);
408 // there were no indices to contract
409 // sort delta8,f,d,T(rl=0),T(rl=1),...,ONE(rl=0),ONE(rl=1),...,unknown
410 // (if there is at least one unknown object, all Ts will be unknown to not change the order)
415 vector<exvector> Tvecs;
416 Tvecs.resize(MAX_REPRESENTATION_LABELS);
417 vector<exvector> ONEvecs;
418 ONEvecs.resize(MAX_REPRESENTATION_LABELS);
421 split_color_string_in_parts(v,delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec);
423 // d_{a,k,l} f_{b,k,l}=0 (includes case a=b)
424 if ((dvec.size()>=1)&&(fvec.size()>=1)) {
425 for (exvector::iterator it1=dvec.begin(); it1!=dvec.end(); ++it1) {
426 for (exvector::iterator it2=fvec.begin(); it2!=fvec.end(); ++it2) {
427 ASSERT(is_ex_exactly_of_type(*it1,color));
428 ASSERT(is_ex_exactly_of_type(*it2,color));
429 color const & col1=ex_to_color(*it1);
430 color const & col2=ex_to_color(*it2);
431 exvector iv_intersect=idx_intersect(col1.seq,col2.seq);
432 if (iv_intersect.size()>=2) return exZERO();
437 // d_{a,k,l} d_{b,k,l}=5/3 delta8_{a,b} (includes case a=b)
438 if (dvec.size()>=2) {
439 for (exvector::iterator it1=dvec.begin(); it1!=dvec.end()-1; ++it1) {
440 for (exvector::iterator it2=it1+1; it2!=dvec.end(); ++it2) {
441 ASSERT(is_ex_exactly_of_type(*it1,color));
442 ASSERT(is_ex_exactly_of_type(*it2,color));
443 color const & col1=ex_to_color(*it1);
444 color const & col2=ex_to_color(*it2);
445 exvector iv_intersect=idx_intersect(col1.seq,col2.seq);
446 if (iv_intersect.size()>=2) {
447 if (iv_intersect.size()==3) {
448 *it1=numeric(40)/numeric(3);
451 int sig1, sig2; // unimportant, since symmetric
452 ex idx1=permute_free_index_to_front(col1.seq,iv_intersect,false,&sig1);
453 ex idx2=permute_free_index_to_front(col2.seq,iv_intersect,false,&sig2);
454 *it1=numeric(5)/numeric(3)*color(color_delta8,idx1,idx2);
457 return nonsimplified_ncmul(recombine_color_string(
458 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
464 // f_{a,k,l} f_{b,k,l}=3 delta8_{a,b} (includes case a=b)
465 if (fvec.size()>=2) {
466 for (exvector::iterator it1=fvec.begin(); it1!=fvec.end()-1; ++it1) {
467 for (exvector::iterator it2=it1+1; it2!=fvec.end(); ++it2) {
468 ASSERT(is_ex_exactly_of_type(*it1,color));
469 ASSERT(is_ex_exactly_of_type(*it2,color));
470 color const & col1=ex_to_color(*it1);
471 color const & col2=ex_to_color(*it2);
472 exvector iv_intersect=idx_intersect(col1.seq,col2.seq);
473 if (iv_intersect.size()>=2) {
474 if (iv_intersect.size()==3) {
479 ex idx1=permute_free_index_to_front(col1.seq,iv_intersect,true,&sig1);
480 ex idx2=permute_free_index_to_front(col2.seq,iv_intersect,true,&sig2);
481 *it1=numeric(sig1*sig2*5)/numeric(3)*color(color_delta8,idx1,idx2);
484 return nonsimplified_ncmul(recombine_color_string(
485 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
491 // d_{a,b,c} T_b T_c = 5/6 T_a
492 // f_{a,b,c} T_b T_c = 3/2 I T_a
493 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
494 if ((Tvecs[rl].size()>=2)&&((dvec.size()>=1)||(fvec.size()>=1))) {
495 for (exvector::iterator it1=Tvecs[rl].begin(); it1!=Tvecs[rl].end()-1; ++it1) {
497 ASSERT(is_ex_exactly_of_type(*it1,color)&&ex_to_color(*it1).type==color_T);
498 ASSERT(is_ex_exactly_of_type(*(it1+1),color)&&ex_to_color(*(it1+1)).type==color_T);
499 iv.push_back(ex_to_color(*it1).seq[0]);
500 iv.push_back(ex_to_color(*(it1+1)).seq[0]);
502 // d_{a,b,c} T_b T_c = 5/6 T_a
503 for (exvector::iterator it2=dvec.begin(); it2!=dvec.end(); ++it2) {
504 ASSERT(is_ex_exactly_of_type(*it2,color)&&ex_to_color(*it2).type==color_d);
505 color const & dref=ex_to_color(*it2);
506 exvector iv_intersect=idx_intersect(dref.seq,iv);
507 if (iv_intersect.size()==2) {
508 int sig; // unimportant, since symmetric
509 ex free_idx=permute_free_index_to_front(dref.seq,iv,false,&sig);
510 *it1=color(color_T,free_idx,rl);
511 *(it1+1)=color(color_ONE,rl);
512 *it2=numeric(5)/numeric(6);
513 return nonsimplified_ncmul(recombine_color_string(
514 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
518 // f_{a,b,c} T_b T_c = 3/2 I T_a
519 for (exvector::iterator it2=fvec.begin(); it2!=fvec.end(); ++it2) {
520 ASSERT(is_ex_exactly_of_type(*it2,color)&&ex_to_color(*it2).type==color_f);
521 color const & fref=ex_to_color(*it2);
522 exvector iv_intersect=idx_intersect(fref.seq,iv);
523 if (iv_intersect.size()==2) {
525 ex free_idx=permute_free_index_to_front(fref.seq,iv,true,&sig);
526 *it1=color(color_T,free_idx,rl);
527 *(it1+1)=color(color_ONE,rl);
528 *it2=numeric(sig*3)/numeric(2)*I;
529 return nonsimplified_ncmul(recombine_color_string(
530 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
537 // clear all ONEs when there is at least one corresponding color_T
538 // in this representation, retain one ONE otherwise
539 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
540 if (Tvecs[rl].size()!=0) {
542 } else if (ONEvecs[rl].size()!=0) {
544 ONEvecs[rl].push_back(color(color_ONE,rl));
548 // return a sorted vector
549 return simplified_ncmul(recombine_color_string(delta8vec,fvec,dvec,Tvecs,
550 ONEvecs,unknownvec));
553 ex color::thisexprseq(exvector const & v) const
555 return color(type,v,representation_label);
558 ex color::thisexprseq(exvector * vp) const
560 return color(type,vp,representation_label);
563 bool color::all_of_type_coloridx(void) const
565 // used only inside of ASSERTs
566 for (exvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
567 if (!is_ex_of_type(*cit,coloridx)) return false;
573 // virtual functions which can be overridden by derived classes
579 // non-virtual functions in this class
583 // static member variables
592 const color some_color;
593 type_info const & typeid_color=typeid(some_color);
599 color color_ONE(unsigned const rl)
601 return color(color::color_ONE,rl);
604 color color_T(ex const & a, unsigned const rl)
606 return color(color::color_T,a,rl);
609 color color_f(ex const & a, ex const & b, ex const & c)
611 return color(color::color_f,a,b,c);
614 color color_d(ex const & a, ex const & b, ex const & c)
616 return color(color::color_d,a,b,c);
619 ex color_h(ex const & a, ex const & b, ex const & c)
621 return color(color::color_d,a,b,c)+I*color(color::color_f,a,b,c);
624 color color_delta8(ex const & a, ex const & b)
626 return color(color::color_delta8,a,b);
629 void split_color_string_in_parts(exvector const & v, exvector & delta8vec,
630 exvector & fvec, exvector & dvec,
631 vector<exvector> & Tvecs,
632 vector<exvector> & ONEvecs,
633 exvector & unknownvec)
635 // if not all elements are of type color, put all Ts in unknownvec to
636 // retain the ordering
638 for (exvector::const_iterator cit=v.begin(); cit!=v.end(); ++cit) {
639 if (!is_ex_exactly_of_type(*cit,color)) {
645 for (exvector::const_iterator cit=v.begin(); cit!=v.end(); ++cit) {
646 if (is_ex_exactly_of_type(*cit,color)) {
647 switch (ex_to_color(*cit).type) {
648 case color::color_delta8:
649 delta8vec.push_back(*cit);
652 fvec.push_back(*cit);
655 dvec.push_back(*cit);
658 ASSERT(ex_to_color(*cit).representation_label<MAX_REPRESENTATION_LABELS);
660 Tvecs[ex_to_color(*cit).representation_label].push_back(*cit);
662 unknownvec.push_back(*cit);
665 case color::color_ONE:
666 ASSERT(ex_to_color(*cit).representation_label<MAX_REPRESENTATION_LABELS);
667 ONEvecs[ex_to_color(*cit).representation_label].push_back(*cit);
670 throw(std::logic_error("invalid type in split_color_string_in_parts()"));
673 unknownvec.push_back(*cit);
678 exvector recombine_color_string(exvector & delta8vec, exvector & fvec,
679 exvector & dvec, vector<exvector> & Tvecs,
680 vector<exvector> & ONEvecs, exvector & unknownvec)
682 unsigned sz=delta8vec.size()+fvec.size()+dvec.size()+unknownvec.size();
683 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
684 sz += Tvecs[rl].size();
685 sz += ONEvecs[rl].size();
690 append_exvector_to_exvector(v,delta8vec);
691 append_exvector_to_exvector(v,fvec);
692 append_exvector_to_exvector(v,dvec);
693 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
694 append_exvector_to_exvector(v,Tvecs[rl]);
695 append_exvector_to_exvector(v,ONEvecs[rl]);
697 append_exvector_to_exvector(v,unknownvec);
701 ex color_trace_of_one_representation_label(exvector const & v)
704 return numeric(COLOR_THREE);
705 } else if (v.size()==1) {
706 ASSERT(is_ex_exactly_of_type(*(v.begin()),color));
710 ex last_element=v1.back();
711 ASSERT(is_ex_exactly_of_type(last_element,color));
712 ASSERT(ex_to_color(last_element).type==color::color_T);
714 ex next_to_last_element=v1.back();
715 ASSERT(is_ex_exactly_of_type(next_to_last_element,color));
716 ASSERT(ex_to_color(next_to_last_element).type==color::color_T);
720 ex const & last_index=ex_to_color(last_element).seq[0];
721 ex const & next_to_last_index=ex_to_color(next_to_last_element).seq[0];
722 ex summation_index=coloridx();
724 v2.push_back(color_T(summation_index)); // don't care about the representation_label
726 // check this formula for SU(N) with N!=3 !!!!!!!!!
727 return numeric(1)/numeric(2*COLOR_THREE)*color_delta8(next_to_last_index,last_index)
728 % color_trace_of_one_representation_label(v1)
729 +numeric(1)/numeric(2)*color_h(next_to_last_index,last_index,summation_index)
730 % color_trace_of_one_representation_label(v2);
732 ex term1=numeric(1)/numeric(2*COLOR_THREE)*color_delta8(next_to_last_index,last_index)
733 % color_trace_of_one_representation_label(v1);
734 cout << "term 1 of trace of " << v.size() << " ts=" << term1 << endl;
735 ex term2=numeric(1)/numeric(2)*color_h(next_to_last_index,last_index,summation_index)
736 % color_trace_of_one_representation_label(v2);
737 cout << "term 2 of trace of " << v.size() << " ts=" << term2 << endl;
742 ex color_trace(exvector const & v, unsigned const rl)
744 ASSERT(rl<MAX_REPRESENTATION_LABELS);
747 v_rest.reserve(v.size()+1); // max size if trace is empty
752 vector<exvector> Tvecs;
753 Tvecs.resize(MAX_REPRESENTATION_LABELS);
754 vector<exvector> ONEvecs;
755 ONEvecs.resize(MAX_REPRESENTATION_LABELS);
758 split_color_string_in_parts(v,delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec);
760 if (unknownvec.size()!=0) {
761 throw(std::invalid_argument("color_trace(): expression must be expanded"));
764 append_exvector_to_exvector(v_rest,delta8vec);
765 append_exvector_to_exvector(v_rest,fvec);
766 append_exvector_to_exvector(v_rest,dvec);
767 for (unsigned i=0; i<MAX_REPRESENTATION_LABELS; ++i) {
769 append_exvector_to_exvector(v_rest,Tvecs[i]);
770 append_exvector_to_exvector(v_rest,ONEvecs[i]);
772 if (Tvecs[i].size()!=0) {
773 v_rest.push_back(color_trace_of_one_representation_label(Tvecs[i]));
774 } else if (ONEvecs[i].size()!=0) {
775 v_rest.push_back(numeric(COLOR_THREE));
777 throw(std::logic_error("color_trace(): representation_label not in color string"));
782 return nonsimplified_ncmul(v_rest);
785 ex simplify_pure_color_string(ex const & e)
787 ASSERT(is_ex_exactly_of_type(e,ncmul));
792 vector<exvector> Tvecs;
793 Tvecs.resize(MAX_REPRESENTATION_LABELS);
794 vector<exvector> ONEvecs;
795 ONEvecs.resize(MAX_REPRESENTATION_LABELS);
798 split_color_string_in_parts(ex_to_ncmul(e).get_factors(),delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec);
800 // search for T_k S T_k (=1/2 Tr(S) - 1/6 S)
801 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
802 if (Tvecs[rl].size()>=2) {
803 for (unsigned i=0; i<Tvecs[rl].size()-1; ++i) {
804 for (unsigned j=i+1; j<Tvecs[rl].size(); ++j) {
805 ex & t1=Tvecs[rl][i];
806 ex & t2=Tvecs[rl][j];
807 ASSERT(is_ex_exactly_of_type(t1,color)&&
808 (ex_to_color(t1).type==color::color_T)&&
809 (ex_to_color(t1).seq.size()==1));
810 ASSERT(is_ex_exactly_of_type(t2,color)&&
811 (ex_to_color(t2).type==color::color_T)&&
812 (ex_to_color(t2).seq.size()==1));
813 coloridx const & idx1=ex_to_coloridx(ex_to_color(t1).seq[0]);
814 coloridx const & idx2=ex_to_coloridx(ex_to_color(t2).seq[0]);
816 if (idx1.is_equal(idx2) && idx1.is_symbolic()) {
818 for (unsigned k=i+1; k<j; ++k) {
819 S.push_back(Tvecs[rl][k]);
823 ex term1=numeric(-1)/numeric(6)*nonsimplified_ncmul(recombine_color_string(
824 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
825 for (unsigned k=i+1; k<j; ++k) {
826 S.push_back(exONE());
828 t1=color_trace_of_one_representation_label(S);
829 ex term2=numeric(1)/numeric(2)*nonsimplified_ncmul(recombine_color_string(
830 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
831 return simplify_color(term1+term2);
838 // TODO: higher contractions!!!!!!!!!!!!!
843 ex simplify_color(ex const & e)
845 // all simplification is done on expanded objects
846 ex e_expanded=e.expand();
848 // simplification of sum=sum of simplifications
849 if (is_ex_exactly_of_type(e_expanded,add)) {
851 for (int i=0; i<e_expanded.nops(); ++i) {
852 sum += simplify_color(e_expanded.op(i));
857 // simplification of commutative product=commutative product of simplifications
858 if (is_ex_exactly_of_type(e_expanded,mul)) {
860 for (int i=0; i<e_expanded.nops(); ++i) {
861 prod *= simplify_color(e_expanded.op(i));
866 // simplification of noncommutative product: test if everything is color
867 if (is_ex_exactly_of_type(e_expanded,ncmul)) {
869 for (int i=0; i<e_expanded.nops(); ++i) {
870 if (!is_ex_exactly_of_type(e_expanded.op(i),color)) {
876 return simplify_pure_color_string(e_expanded);
880 // cannot do anything
884 ex brute_force_sum_color_indices(ex const & e)
886 exvector iv_all=e.get_indices();
889 // find double symbolic indices
890 if (iv_all.size()<2) return e;
891 for (exvector::const_iterator cit1=iv_all.begin(); cit1!=iv_all.end()-1; ++cit1) {
892 ASSERT(is_ex_of_type(*cit1,coloridx));
893 for (exvector::const_iterator cit2=cit1+1; cit2!=iv_all.end(); ++cit2) {
894 ASSERT(is_ex_of_type(*cit2,coloridx));
895 if (ex_to_coloridx(*cit1).is_symbolic() &&
896 ex_to_coloridx(*cit1).is_equal(ex_to_coloridx(*cit2))) {
897 iv_double.push_back(*cit1);
904 counter.resize(iv_double.size());
906 for (l=0; unsigned(l)<iv_double.size(); ++l) {
914 for (l=0; unsigned(l)<iv_double.size(); ++l) {
915 term=term.subs(iv_double[l]==coloridx((unsigned)(counter[l])));
916 //iv_double[l].print(cout);
917 //cout << " " << counter[l] << " ";
922 // increment counter[]
923 l=iv_double.size()-1;
924 while ((l>=0)&&((++counter[l])>COLOR_EIGHT)) {
928 if (l<2) { cout << counter[0] << counter[1] << endl; }
935 void append_exvector_to_exvector(exvector & dest, exvector const & source)
937 for (exvector::const_iterator cit=source.begin(); cit!=source.end(); ++cit) {
938 dest.push_back(*cit);