3 * Implementation of GiNaC's color objects.
4 * No real implementation yet, to be done. */
7 * GiNaC Copyright (C) 1999-2000 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"
40 #ifndef NO_NAMESPACE_GINAC
42 #endif // ndef NO_NAMESPACE_GINAC
44 GINAC_IMPLEMENT_REGISTERED_CLASS(color, indexed)
47 // default constructor, destructor, copy constructor assignment operator and helpers
52 color::color() : type(invalid), representation_label(0)
54 debugmsg("color default constructor",LOGLEVEL_CONSTRUCT);
55 tinfo_key=TINFO_color;
60 debugmsg("color destructor",LOGLEVEL_DESTRUCT);
64 color::color(const color & other)
66 debugmsg("color copy constructor",LOGLEVEL_CONSTRUCT);
70 const color & color::operator=(const color & other)
72 debugmsg("color operator=",LOGLEVEL_ASSIGNMENT);
82 void color::copy(const color & other)
84 inherited::copy(other);
86 representation_label=other.representation_label;
89 void color::destroy(bool call_parent)
92 inherited::destroy(call_parent);
102 color::color(color_types const t, unsigned rl) : type(t), representation_label(rl)
104 debugmsg("color constructor from color_types,unsigned",LOGLEVEL_CONSTRUCT);
105 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
106 tinfo_key=TINFO_color;
107 GINAC_ASSERT(all_of_type_coloridx());
110 color::color(color_types const t, const ex & i1, unsigned rl)
111 : inherited(i1), type(t), representation_label(rl)
113 debugmsg("color constructor from color_types,ex,unsigned",LOGLEVEL_CONSTRUCT);
114 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
115 tinfo_key=TINFO_color;
116 GINAC_ASSERT(all_of_type_coloridx());
119 color::color(color_types const t, const ex & i1, const ex & i2, unsigned rl)
120 : inherited(i1,i2), type(t), representation_label(rl)
122 debugmsg("color constructor from color_types,ex,ex,unsigned",LOGLEVEL_CONSTRUCT);
123 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
124 tinfo_key=TINFO_color;
125 GINAC_ASSERT(all_of_type_coloridx());
128 color::color(color_types const t, const ex & i1, const ex & i2, const ex & i3,
129 unsigned rl) : inherited(i1,i2,i3), type(t), representation_label(rl)
131 debugmsg("color constructor from color_types,ex,ex,ex,unsigned",LOGLEVEL_CONSTRUCT);
132 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
133 tinfo_key=TINFO_color;
134 GINAC_ASSERT(all_of_type_coloridx());
137 color::color(color_types const t, const exvector & iv, unsigned rl)
138 : inherited(iv), type(t), representation_label(rl)
140 debugmsg("color constructor from color_types,exvector,unsigned",LOGLEVEL_CONSTRUCT);
141 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
142 tinfo_key=TINFO_color;
143 GINAC_ASSERT(all_of_type_coloridx());
146 color::color(color_types const t, exvector * ivp, unsigned rl)
147 : inherited(ivp), type(t), representation_label(rl)
149 debugmsg("color constructor from color_types,exvector *,unsigned",LOGLEVEL_CONSTRUCT);
150 GINAC_ASSERT(representation_label<MAX_REPRESENTATION_LABELS);
151 tinfo_key=TINFO_color;
152 GINAC_ASSERT(all_of_type_coloridx());
159 /** Construct object from archive_node. */
160 color::color(const archive_node &n, const lst &sym_lst) : inherited(n, sym_lst)
162 debugmsg("color constructor from archive_node", LOGLEVEL_CONSTRUCT);
164 if (!(n.find_unsigned("type", ty)))
165 throw (std::runtime_error("unknown color type in archive"));
166 type = (color_types)ty;
167 if (!(n.find_unsigned("representation", representation_label)))
168 throw (std::runtime_error("unknown color representation label in archive"));
171 /** Unarchive the object. */
172 ex color::unarchive(const archive_node &n, const lst &sym_lst)
174 return (new color(n, sym_lst))->setflag(status_flags::dynallocated);
177 /** Archive the object. */
178 void color::archive(archive_node &n) const
180 inherited::archive(n);
181 n.add_unsigned("type", type);
182 n.add_unsigned("representation", representation_label);
186 // functions overriding virtual functions from bases classes
191 basic * color::duplicate() const
193 debugmsg("color duplicate",LOGLEVEL_DUPLICATE);
194 return new color(*this);
197 void color::printraw(ostream & os) const
199 debugmsg("color printraw",LOGLEVEL_PRINT);
200 os << "color(type=" << (unsigned)type
201 << ",representation_label=" << representation_label
204 os << ",hash=" << hashvalue << ",flags=" << flags << ")";
207 void color::printtree(ostream & os, unsigned indent) const
209 debugmsg("color printtree",LOGLEVEL_PRINT);
210 os << string(indent,' ') << "color object: "
211 << "type=" << (unsigned)type
212 << ",representation_label=" << representation_label << ", ";
213 os << seq.size() << " indices" << endl;
214 printtreeindices(os,indent);
215 os << string(indent,' ') << "hash=" << hashvalue
216 << " (0x" << hex << hashvalue << dec << ")"
217 << ", flags=" << flags << endl;
220 void color::print(ostream & os, unsigned upper_precedence) const
222 debugmsg("color print",LOGLEVEL_PRINT);
226 if (representation_label!=0) {
227 os << "^(" << representation_label << ")";
244 os << "INVALID_COLOR_OBJECT";
250 void color::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const
252 debugmsg("color print csrc",LOGLEVEL_PRINT);
253 print(os,upper_precedence);
256 bool color::info(unsigned inf) const
258 return inherited::info(inf);
261 #define CMPINDICES(A,B,C) ((idx1.get_value()==(A))&&(idx2.get_value()==(B))&&(idx3.get_value()==(C)))
263 ex color::eval(int level) const
265 // canonicalize indices
267 bool antisymmetric=false;
271 antisymmetric=true; // no break here!
276 int sig=canonicalize_indices(iv,antisymmetric);
278 // something has changed while sorting indices, more evaluations later
279 if (sig==0) return _ex0();
280 return ex(sig)*color(type,iv,representation_label);
285 // nothing to canonicalize
292 GINAC_ASSERT(seq.size()==2);
293 const coloridx & idx1=ex_to_coloridx(seq[0]);
294 const coloridx & idx2=ex_to_coloridx(seq[1]);
296 // check for delta8_{a,a} where a is a symbolic index, replace by 8
297 if ((idx1.is_symbolic())&&(idx1.is_equal_same_type(idx2))) {
298 return ex(COLOR_EIGHT);
301 // check for delta8_{a,b} where a and b are numeric indices, replace by 0 or 1
302 if ((!idx1.is_symbolic())&&(!idx2.is_symbolic())) {
303 if ((idx1.get_value()!=idx2.get_value())) {
312 // check for d_{a,a,c} (=0) when a is symbolic
314 GINAC_ASSERT(seq.size()==3);
315 const coloridx & idx1=ex_to_coloridx(seq[0]);
316 const coloridx & idx2=ex_to_coloridx(seq[1]);
317 const coloridx & idx3=ex_to_coloridx(seq[2]);
319 if (idx1.is_equal_same_type(idx2) && idx1.is_symbolic()) {
321 } else if (idx2.is_equal_same_type(idx3) && idx2.is_symbolic()) {
325 // check for three numeric indices
326 if (!(idx1.is_symbolic()||idx2.is_symbolic()||idx3.is_symbolic())) {
327 GINAC_ASSERT(idx1.get_value()<=idx2.get_value());
328 GINAC_ASSERT(idx2.get_value()<=idx3.get_value());
329 if (CMPINDICES(1,4,6)||CMPINDICES(1,5,7)||CMPINDICES(2,5,6)||
330 CMPINDICES(3,4,4)||CMPINDICES(3,5,5)) {
332 } else if (CMPINDICES(2,4,7)||CMPINDICES(3,6,6)||CMPINDICES(3,7,7)) {
334 } else if (CMPINDICES(1,1,8)||CMPINDICES(2,2,8)||CMPINDICES(3,3,8)) {
335 return 1/sqrt(numeric(3));
336 } else if (CMPINDICES(8,8,8)) {
337 return -1/sqrt(numeric(3));
338 } else if (CMPINDICES(4,4,8)||CMPINDICES(5,5,8)||CMPINDICES(6,6,8)||CMPINDICES(7,7,8)) {
339 return -1/(2*sqrt(numeric(3)));
347 GINAC_ASSERT(seq.size()==3);
348 const coloridx & idx1=ex_to_coloridx(seq[0]);
349 const coloridx & idx2=ex_to_coloridx(seq[1]);
350 const coloridx & idx3=ex_to_coloridx(seq[2]);
352 // check for three numeric indices
353 if (!(idx1.is_symbolic()||idx2.is_symbolic()||idx3.is_symbolic())) {
354 GINAC_ASSERT(idx1.get_value()<=idx2.get_value());
355 GINAC_ASSERT(idx2.get_value()<=idx3.get_value());
356 if (CMPINDICES(1,2,3)) {
358 } else if (CMPINDICES(1,4,7)||CMPINDICES(2,4,6)||
359 CMPINDICES(2,5,7)||CMPINDICES(3,4,5)) {
361 } else if (CMPINDICES(1,5,6)||CMPINDICES(3,6,7)) {
363 } else if (CMPINDICES(4,5,8)||CMPINDICES(6,7,8)) {
364 return sqrt(numeric(3))/2;
365 } else if (CMPINDICES(8,8,8)) {
366 return -1/sqrt(numeric(3));
367 } else if (CMPINDICES(4,4,8)||CMPINDICES(5,5,8)||CMPINDICES(6,6,8)||CMPINDICES(7,7,8)) {
368 return -1/(2*sqrt(numeric(3)));
375 // nothing to evaluate
384 int color::compare_same_type(const basic & other) const
386 GINAC_ASSERT(other.tinfo() == TINFO_color);
387 const color *o = static_cast<const color *>(&other);
389 if (representation_label==o->representation_label) {
390 return inherited::compare_same_type(other);
392 return representation_label < o->representation_label ? -1 : 1;
394 return type < o->type ? -1 : 1;
397 bool color::is_equal_same_type(const basic & other) const
399 GINAC_ASSERT(other.tinfo() == TINFO_color);
400 const color *o = static_cast<const color *>(&other);
401 if (type!=o->type) return false;
402 if (representation_label!=o->representation_label) return false;
403 return inherited::is_equal_same_type(other);
408 ex color::simplify_ncmul(const exvector & v) const
410 // simplifications: contract delta8_{a,b} where possible
411 // sort delta8,f,d,T(rl=0),T(rl=1),...,ONE(rl=0),ONE(rl=1),...
412 // remove superfluous ONEs
414 // contract indices of delta8_{a,b} if they are different and symbolic
416 exvector v_contracted=v;
417 unsigned replacements;
418 bool something_changed=false;
420 exvector::iterator it=v_contracted.begin();
421 while (it!=v_contracted.end()) {
422 // process only delta8 objects
423 if (is_ex_exactly_of_type(*it,color) && (ex_to_color(*it).type==color_delta8)) {
424 color & d8=ex_to_nonconst_color(*it);
425 GINAC_ASSERT(d8.seq.size()==2);
426 const coloridx & first_idx=ex_to_coloridx(d8.seq[0]);
427 const coloridx & second_idx=ex_to_coloridx(d8.seq[1]);
428 // delta8_{a,a} should have been contracted in color::eval()
429 GINAC_ASSERT((!first_idx.is_equal(second_idx))||(!first_idx.is_symbolic()));
430 ex saved_delta8=*it; // save to restore it later
432 // try to contract first index
434 if (first_idx.is_symbolic()) {
435 replacements = subs_index_in_exvector(v_contracted,first_idx,second_idx);
436 if (replacements==1) {
437 // not contracted except in itself, restore delta8 object
440 // a contracted index should occur exactly twice
441 GINAC_ASSERT(replacements==2);
443 something_changed=true;
447 // try second index only if first was not contracted
448 if ((replacements==1)&&(second_idx.is_symbolic())) {
449 // first index not contracted, *it is guaranteed to be the original delta8 object
450 replacements = subs_index_in_exvector(v_contracted,second_idx,first_idx);
451 if (replacements==1) {
452 // not contracted except in itself, restore delta8 object
455 // a contracted index should occur exactly twice
456 GINAC_ASSERT(replacements==2);
458 something_changed=true;
465 if (something_changed) {
466 // do more simplifications later
467 return nonsimplified_ncmul(v_contracted);
470 // there were no indices to contract
471 // sort delta8,f,d,T(rl=0),T(rl=1),...,ONE(rl=0),ONE(rl=1),...,unknown
472 // (if there is at least one unknown object, all Ts will be unknown to not change the order)
477 exvectorvector Tvecs;
478 Tvecs.resize(MAX_REPRESENTATION_LABELS);
479 exvectorvector ONEvecs;
480 ONEvecs.resize(MAX_REPRESENTATION_LABELS);
483 split_color_string_in_parts(v,delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec);
485 // d_{a,k,l} f_{b,k,l}=0 (includes case a=b)
486 if ((dvec.size()>=1)&&(fvec.size()>=1)) {
487 for (exvector::iterator it1=dvec.begin(); it1!=dvec.end(); ++it1) {
488 for (exvector::iterator it2=fvec.begin(); it2!=fvec.end(); ++it2) {
489 GINAC_ASSERT(is_ex_exactly_of_type(*it1,color));
490 GINAC_ASSERT(is_ex_exactly_of_type(*it2,color));
491 const color & col1=ex_to_color(*it1);
492 const color & col2=ex_to_color(*it2);
493 exvector iv_intersect=idx_intersect(col1.seq,col2.seq);
494 if (iv_intersect.size()>=2) return _ex0();
499 // d_{a,k,l} d_{b,k,l}=5/3 delta8_{a,b} (includes case a=b)
500 if (dvec.size()>=2) {
501 for (exvector::iterator it1=dvec.begin(); it1!=dvec.end()-1; ++it1) {
502 for (exvector::iterator it2=it1+1; it2!=dvec.end(); ++it2) {
503 GINAC_ASSERT(is_ex_exactly_of_type(*it1,color));
504 GINAC_ASSERT(is_ex_exactly_of_type(*it2,color));
505 const color & col1=ex_to_color(*it1);
506 const color & col2=ex_to_color(*it2);
507 exvector iv_intersect=idx_intersect(col1.seq,col2.seq);
508 if (iv_intersect.size()>=2) {
509 if (iv_intersect.size()==3) {
510 *it1=numeric(40)/numeric(3);
513 int sig1, sig2; // unimportant, since symmetric
514 ex idx1=permute_free_index_to_front(col1.seq,iv_intersect,false,&sig1);
515 ex idx2=permute_free_index_to_front(col2.seq,iv_intersect,false,&sig2);
516 *it1=numeric(5)/numeric(3)*color(color_delta8,idx1,idx2);
519 return nonsimplified_ncmul(recombine_color_string(
520 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
526 // f_{a,k,l} f_{b,k,l}=3 delta8_{a,b} (includes case a=b)
527 if (fvec.size()>=2) {
528 for (exvector::iterator it1=fvec.begin(); it1!=fvec.end()-1; ++it1) {
529 for (exvector::iterator it2=it1+1; it2!=fvec.end(); ++it2) {
530 GINAC_ASSERT(is_ex_exactly_of_type(*it1,color));
531 GINAC_ASSERT(is_ex_exactly_of_type(*it2,color));
532 const color & col1=ex_to_color(*it1);
533 const color & col2=ex_to_color(*it2);
534 exvector iv_intersect=idx_intersect(col1.seq,col2.seq);
535 if (iv_intersect.size()>=2) {
536 if (iv_intersect.size()==3) {
541 ex idx1=permute_free_index_to_front(col1.seq,iv_intersect,true,&sig1);
542 ex idx2=permute_free_index_to_front(col2.seq,iv_intersect,true,&sig2);
543 *it1=numeric(sig1*sig2*5)/numeric(3)*color(color_delta8,idx1,idx2);
546 return nonsimplified_ncmul(recombine_color_string(
547 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
553 // d_{a,b,c} T_b T_c = 5/6 T_a
554 // f_{a,b,c} T_b T_c = 3/2 I T_a
555 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
556 if ((Tvecs[rl].size()>=2)&&((dvec.size()>=1)||(fvec.size()>=1))) {
557 for (exvector::iterator it1=Tvecs[rl].begin(); it1!=Tvecs[rl].end()-1; ++it1) {
559 GINAC_ASSERT(is_ex_exactly_of_type(*it1,color)&&ex_to_color(*it1).type==color_T);
560 GINAC_ASSERT(is_ex_exactly_of_type(*(it1+1),color)&&ex_to_color(*(it1+1)).type==color_T);
561 iv.push_back(ex_to_color(*it1).seq[0]);
562 iv.push_back(ex_to_color(*(it1+1)).seq[0]);
564 // d_{a,b,c} T_b T_c = 5/6 T_a
565 for (exvector::iterator it2=dvec.begin(); it2!=dvec.end(); ++it2) {
566 GINAC_ASSERT(is_ex_exactly_of_type(*it2,color)&&ex_to_color(*it2).type==color_d);
567 const color & dref=ex_to_color(*it2);
568 exvector iv_intersect=idx_intersect(dref.seq,iv);
569 if (iv_intersect.size()==2) {
570 int sig; // unimportant, since symmetric
571 ex free_idx=permute_free_index_to_front(dref.seq,iv,false,&sig);
572 *it1=color(color_T,free_idx,rl);
573 *(it1+1)=color(color_ONE,rl);
574 *it2=numeric(5)/numeric(6);
575 return nonsimplified_ncmul(recombine_color_string(
576 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
580 // f_{a,b,c} T_b T_c = 3/2 I T_a
581 for (exvector::iterator it2=fvec.begin(); it2!=fvec.end(); ++it2) {
582 GINAC_ASSERT(is_ex_exactly_of_type(*it2,color)&&ex_to_color(*it2).type==color_f);
583 const color & fref=ex_to_color(*it2);
584 exvector iv_intersect=idx_intersect(fref.seq,iv);
585 if (iv_intersect.size()==2) {
587 ex free_idx=permute_free_index_to_front(fref.seq,iv,true,&sig);
588 *it1=color(color_T,free_idx,rl);
589 *(it1+1)=color(color_ONE,rl);
590 *it2=numeric(sig*3)/numeric(2)*I;
591 return nonsimplified_ncmul(recombine_color_string(
592 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
599 // clear all ONEs when there is at least one corresponding color_T
600 // in this representation, retain one ONE otherwise
601 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
602 if (Tvecs[rl].size()!=0) {
604 } else if (ONEvecs[rl].size()!=0) {
606 ONEvecs[rl].push_back(color(color_ONE,rl));
610 // return a sorted vector
611 return simplified_ncmul(recombine_color_string(delta8vec,fvec,dvec,Tvecs,
612 ONEvecs,unknownvec));
615 ex color::thisexprseq(const exvector & v) const
617 return color(type,v,representation_label);
620 ex color::thisexprseq(exvector * vp) const
622 return color(type,vp,representation_label);
625 bool color::all_of_type_coloridx(void) const
627 // used only inside of ASSERTs
628 for (exvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) {
629 if (!is_ex_of_type(*cit,coloridx)) return false;
635 // virtual functions which can be overridden by derived classes
641 // non-virtual functions in this class
645 // static member variables
654 const color some_color;
655 const type_info & typeid_color=typeid(some_color);
661 color color_ONE(unsigned rl)
663 return color(color::color_ONE,rl);
666 color color_T(const ex & a, unsigned rl)
668 return color(color::color_T,a,rl);
671 color color_f(const ex & a, const ex & b, const ex & c)
673 return color(color::color_f,a,b,c);
676 color color_d(const ex & a, const ex & b, const ex & c)
678 return color(color::color_d,a,b,c);
681 ex color_h(const ex & a, const ex & b, const ex & c)
683 return color(color::color_d,a,b,c)+I*color(color::color_f,a,b,c);
686 color color_delta8(const ex & a, const ex & b)
688 return color(color::color_delta8,a,b);
691 void split_color_string_in_parts(const exvector & v, exvector & delta8vec,
692 exvector & fvec, exvector & dvec,
693 exvectorvector & Tvecs,
694 exvectorvector & ONEvecs,
695 exvector & unknownvec)
697 // if not all elements are of type color, put all Ts in unknownvec to
698 // retain the ordering
700 for (exvector::const_iterator cit=v.begin(); cit!=v.end(); ++cit) {
701 if (!is_ex_exactly_of_type(*cit,color)) {
707 for (exvector::const_iterator cit=v.begin(); cit!=v.end(); ++cit) {
708 if (is_ex_exactly_of_type(*cit,color)) {
709 switch (ex_to_color(*cit).type) {
710 case color::color_delta8:
711 delta8vec.push_back(*cit);
714 fvec.push_back(*cit);
717 dvec.push_back(*cit);
720 GINAC_ASSERT(ex_to_color(*cit).representation_label<MAX_REPRESENTATION_LABELS);
722 Tvecs[ex_to_color(*cit).representation_label].push_back(*cit);
724 unknownvec.push_back(*cit);
727 case color::color_ONE:
728 GINAC_ASSERT(ex_to_color(*cit).representation_label<MAX_REPRESENTATION_LABELS);
729 ONEvecs[ex_to_color(*cit).representation_label].push_back(*cit);
732 throw(std::logic_error("invalid type in split_color_string_in_parts()"));
735 unknownvec.push_back(*cit);
740 exvector recombine_color_string(exvector & delta8vec, exvector & fvec,
741 exvector & dvec, exvectorvector & Tvecs,
742 exvectorvector & ONEvecs, exvector & unknownvec)
744 unsigned sz=delta8vec.size()+fvec.size()+dvec.size()+unknownvec.size();
745 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
746 sz += Tvecs[rl].size();
747 sz += ONEvecs[rl].size();
752 append_exvector_to_exvector(v,delta8vec);
753 append_exvector_to_exvector(v,fvec);
754 append_exvector_to_exvector(v,dvec);
755 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
756 append_exvector_to_exvector(v,Tvecs[rl]);
757 append_exvector_to_exvector(v,ONEvecs[rl]);
759 append_exvector_to_exvector(v,unknownvec);
763 ex color_trace_of_one_representation_label(const exvector & v)
766 return numeric(COLOR_THREE);
767 } else if (v.size()==1) {
768 GINAC_ASSERT(is_ex_exactly_of_type(*(v.begin()),color));
772 ex last_element=v1.back();
773 GINAC_ASSERT(is_ex_exactly_of_type(last_element,color));
774 GINAC_ASSERT(ex_to_color(last_element).type==color::color_T);
776 ex next_to_last_element=v1.back();
777 GINAC_ASSERT(is_ex_exactly_of_type(next_to_last_element,color));
778 GINAC_ASSERT(ex_to_color(next_to_last_element).type==color::color_T);
782 const ex & last_index=ex_to_color(last_element).seq[0];
783 const ex & next_to_last_index=ex_to_color(next_to_last_element).seq[0];
784 ex summation_index=coloridx();
786 v2.push_back(color_T(summation_index)); // don't care about the representation_label
788 // FIXME: check this formula for SU(N) with N!=3
789 return numeric(1)/numeric(2*COLOR_THREE)*color_delta8(next_to_last_index,last_index)
790 % color_trace_of_one_representation_label(v1)
791 +numeric(1)/numeric(2)*color_h(next_to_last_index,last_index,summation_index)
792 % color_trace_of_one_representation_label(v2);
794 ex term1=numeric(1)/numeric(2*COLOR_THREE)*color_delta8(next_to_last_index,last_index)
795 % color_trace_of_one_representation_label(v1);
796 cout << "term 1 of trace of " << v.size() << " ts=" << term1 << endl;
797 ex term2=numeric(1)/numeric(2)*color_h(next_to_last_index,last_index,summation_index)
798 % color_trace_of_one_representation_label(v2);
799 cout << "term 2 of trace of " << v.size() << " ts=" << term2 << endl;
804 ex color_trace(const exvector & v, unsigned rl)
806 GINAC_ASSERT(rl<MAX_REPRESENTATION_LABELS);
809 v_rest.reserve(v.size()+1); // max size if trace is empty
814 exvectorvector Tvecs;
815 Tvecs.resize(MAX_REPRESENTATION_LABELS);
816 exvectorvector ONEvecs;
817 ONEvecs.resize(MAX_REPRESENTATION_LABELS);
820 split_color_string_in_parts(v,delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec);
822 if (unknownvec.size()!=0) {
823 throw(std::invalid_argument("color_trace(): expression must be expanded"));
826 append_exvector_to_exvector(v_rest,delta8vec);
827 append_exvector_to_exvector(v_rest,fvec);
828 append_exvector_to_exvector(v_rest,dvec);
829 for (unsigned i=0; i<MAX_REPRESENTATION_LABELS; ++i) {
831 append_exvector_to_exvector(v_rest,Tvecs[i]);
832 append_exvector_to_exvector(v_rest,ONEvecs[i]);
834 if (Tvecs[i].size()!=0) {
835 v_rest.push_back(color_trace_of_one_representation_label(Tvecs[i]));
836 } else if (ONEvecs[i].size()!=0) {
837 v_rest.push_back(numeric(COLOR_THREE));
839 throw(std::logic_error("color_trace(): representation_label not in color string"));
844 return nonsimplified_ncmul(v_rest);
847 ex simplify_pure_color_string(const ex & e)
849 GINAC_ASSERT(is_ex_exactly_of_type(e,ncmul));
854 exvectorvector Tvecs;
855 Tvecs.resize(MAX_REPRESENTATION_LABELS);
856 exvectorvector ONEvecs;
857 ONEvecs.resize(MAX_REPRESENTATION_LABELS);
860 split_color_string_in_parts(ex_to_ncmul(e).get_factors(),delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec);
862 // search for T_k S T_k (=1/2 Tr(S) - 1/6 S)
863 for (unsigned rl=0; rl<MAX_REPRESENTATION_LABELS; ++rl) {
864 if (Tvecs[rl].size()>=2) {
865 for (unsigned i=0; i<Tvecs[rl].size()-1; ++i) {
866 for (unsigned j=i+1; j<Tvecs[rl].size(); ++j) {
867 ex & t1=Tvecs[rl][i];
868 ex & t2=Tvecs[rl][j];
869 GINAC_ASSERT(is_ex_exactly_of_type(t1,color)&&
870 (ex_to_color(t1).type==color::color_T)&&
871 (ex_to_color(t1).seq.size()==1));
872 GINAC_ASSERT(is_ex_exactly_of_type(t2,color)&&
873 (ex_to_color(t2).type==color::color_T)&&
874 (ex_to_color(t2).seq.size()==1));
875 const coloridx & idx1=ex_to_coloridx(ex_to_color(t1).seq[0]);
876 const coloridx & idx2=ex_to_coloridx(ex_to_color(t2).seq[0]);
878 if (idx1.is_equal(idx2) && idx1.is_symbolic()) {
880 for (unsigned k=i+1; k<j; ++k) {
881 S.push_back(Tvecs[rl][k]);
885 ex term1=numeric(-1)/numeric(6)*nonsimplified_ncmul(recombine_color_string(
886 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
887 for (unsigned k=i+1; k<j; ++k) {
890 t1=color_trace_of_one_representation_label(S);
891 ex term2=numeric(1)/numeric(2)*nonsimplified_ncmul(recombine_color_string(
892 delta8vec,fvec,dvec,Tvecs,ONEvecs,unknownvec));
893 return simplify_color(term1+term2);
900 // FIXME: higher contractions
905 ex simplify_color(const ex & e)
907 // all simplification is done on expanded objects
908 ex e_expanded=e.expand();
910 // simplification of sum=sum of simplifications
911 if (is_ex_exactly_of_type(e_expanded,add)) {
913 for (unsigned i=0; i<e_expanded.nops(); ++i)
914 sum += simplify_color(e_expanded.op(i));
919 // simplification of commutative product=commutative product of simplifications
920 if (is_ex_exactly_of_type(e_expanded,mul)) {
922 for (unsigned i=0; i<e_expanded.nops(); ++i)
923 prod *= simplify_color(e_expanded.op(i));
928 // simplification of noncommutative product: test if everything is color
929 if (is_ex_exactly_of_type(e_expanded,ncmul)) {
931 for (unsigned i=0; i<e_expanded.nops(); ++i) {
932 if (!is_ex_exactly_of_type(e_expanded.op(i),color)) {
938 return simplify_pure_color_string(e_expanded);
942 // cannot do anything
946 ex brute_force_sum_color_indices(const ex & e)
948 exvector iv_all=e.get_indices();
951 // find double symbolic indices
952 if (iv_all.size()<2) return e;
953 for (exvector::const_iterator cit1=iv_all.begin(); cit1!=iv_all.end()-1; ++cit1) {
954 GINAC_ASSERT(is_ex_of_type(*cit1,coloridx));
955 for (exvector::const_iterator cit2=cit1+1; cit2!=iv_all.end(); ++cit2) {
956 GINAC_ASSERT(is_ex_of_type(*cit2,coloridx));
957 if (ex_to_coloridx(*cit1).is_symbolic() &&
958 ex_to_coloridx(*cit1).is_equal(ex_to_coloridx(*cit2))) {
959 iv_double.push_back(*cit1);
966 counter.resize(iv_double.size());
968 for (l=0; unsigned(l)<iv_double.size(); ++l) {
976 for (l=0; unsigned(l)<iv_double.size(); ++l) {
977 term=term.subs(iv_double[l]==coloridx((unsigned)(counter[l])));
978 //iv_double[l].print(cout);
979 //cout << " " << counter[l] << " ";
984 // increment counter[]
985 l=iv_double.size()-1;
986 while ((l>=0)&&((++counter[l])>(int)COLOR_EIGHT)) {
990 if (l<2) { cout << counter[0] << counter[1] << endl; }
997 void append_exvector_to_exvector(exvector & dest, const exvector & source)
999 for (exvector::const_iterator cit=source.begin(); cit!=source.end(); ++cit) {
1000 dest.push_back(*cit);
1004 #ifndef NO_NAMESPACE_GINAC
1005 } // namespace GiNaC
1006 #endif // ndef NO_NAMESPACE_GINAC