ginac/idx.cpp

   1 /** @file idx.cpp
   2  *
   3  *  Implementation of GiNaC's indices. */
   4
   5 /*
   6  *  GiNaC Copyright (C) 1999-2000 Johannes Gutenberg University Mainz, Germany
   7  *
   8  *  This program is free software; you can redistribute it and/or modify
   9  *  it under the terms of the GNU General Public License as published by
  10  *  the Free Software Foundation; either version 2 of the License, or
  11  *  (at your option) any later version.
  12  *
  13  *  This program is distributed in the hope that it will be useful,
  14  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
  15  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  16  *  GNU General Public License for more details.
  17  *
  18  *  You should have received a copy of the GNU General Public License
  19  *  along with this program; if not, write to the Free Software
  20  *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  21  */
  22
  23 #include <stdexcept>
  24
  25 #include "idx.h"
  26 #include "ex.h"
  27 #include "lst.h"
  28 #include "relational.h"
  29 #include "utils.h"
  30 #include "debugmsg.h"
  31
  32 #ifndef NO_GINAC_NAMESPACE
  33 namespace GiNaC {
  34 #endif // ndef NO_GINAC_NAMESPACE
  35
  36 //////////
  37 // default constructor, destructor, copy constructor assignment operator and helpers
  38 //////////
  39
  40 // public
  41
  42 idx::idx() : basic(TINFO_idx), symbolic(true), covariant(false)
  43 {
  44     debugmsg("idx default constructor",LOGLEVEL_CONSTRUCT);
  45     serial=next_serial++;
  46     name="index"+ToString(serial);
  47 }
  48
  49 idx::~idx()
  50 {
  51     debugmsg("idx destructor",LOGLEVEL_DESTRUCT);
  52     destroy(0);
  53 }
  54
  55 idx::idx(const idx & other)
  56 {
  57     debugmsg("idx copy constructor",LOGLEVEL_CONSTRUCT);
  58     copy(other);
  59 }
  60
  61 const idx & idx::operator=(const idx & other)
  62 {
  63     debugmsg("idx operator=",LOGLEVEL_ASSIGNMENT);
  64     if (this != &other) {
  65         destroy(1);
  66         copy(other);
  67     }
  68     return *this;
  69 }
  70
  71 // protected
  72
  73 void idx::copy(const idx & other)
  74 {
  75     basic::copy(other);
  76     serial=other.serial;
  77     symbolic=other.symbolic;
  78     name=other.name;
  79     value=other.value;
  80     covariant=other.covariant;
  81 }
  82
  83 void idx::destroy(bool call_parent)
  84 {
  85     if (call_parent) basic::destroy(call_parent);
  86 }
  87
  88 //////////
  89 // other constructors
  90 //////////
  91
  92 // public
  93
  94 idx::idx(bool cov) : basic(TINFO_idx), symbolic(true), covariant(cov)
  95 {
  96     debugmsg("idx constructor from bool",LOGLEVEL_CONSTRUCT);
  97     serial=next_serial++;
  98     name="index"+ToString(serial);
  99 }
 100
 101 idx::idx(const string & n, bool cov) : basic(TINFO_idx),
 102     symbolic(true), name(n), covariant(cov)
 103 {
 104     debugmsg("idx constructor from string,bool",LOGLEVEL_CONSTRUCT);
 105     serial=next_serial++;
 106 }
 107
 108 idx::idx(const char * n, bool cov) : basic(TINFO_idx),
 109     symbolic(true), name(n), covariant(cov)
 110 {
 111     debugmsg("idx constructor from char*,bool",LOGLEVEL_CONSTRUCT);
 112     serial=next_serial++;
 113 }
 114
 115 idx::idx(unsigned v, bool cov) : basic(TINFO_idx),
 116     symbolic(false), value(v), covariant(cov)
 117 {
 118     debugmsg("idx constructor from unsigned,bool",LOGLEVEL_CONSTRUCT);
 119     serial=0;
 120 }
 121
 122
 123 //////////
 124 // functions overriding virtual functions from bases classes
 125 //////////
 126
 127 // public
 128
 129 basic * idx::duplicate() const
 130 {
 131     debugmsg("idx duplicate",LOGLEVEL_DUPLICATE);
 132     return new idx(*this);
 133 }
 134
 135 void idx::printraw(ostream & os) const
 136 {
 137     debugmsg("idx printraw",LOGLEVEL_PRINT);
 138
 139     os << "idx(";
 140
 141     if (symbolic) {
 142         os << "symbolic,name=" << name;
 143     } else {
 144         os << "non symbolic,value=" << value;
 145     }
 146
 147     if (covariant) {
 148         os << ",covariant";
 149     } else {
 150         os << ",contravariant";
 151     }
 152
 153     os << ",serial=" << serial;
 154     os << ",hash=" << hashvalue << ",flags=" << flags;
 155     os << ")";
 156 }
 157
 158 void idx::printtree(ostream & os, unsigned indent) const
 159 {
 160     debugmsg("idx printtree",LOGLEVEL_PRINT);
 161
 162     os << string(indent,' ') << "idx: ";
 163
 164     if (symbolic) {
 165         os << "symbolic,name=" << name;
 166     } else {
 167         os << "non symbolic,value=" << value;
 168     }
 169
 170     if (covariant) {
 171         os << ",covariant";
 172     } else {
 173         os << ",contravariant";
 174     }
 175
 176     os << ", serial=" << serial
 177        << ", hash=" << hashvalue << " (0x" << hex << hashvalue << dec << ")"
 178        << ", flags=" << flags << endl;
 179 }
 180
 181 void idx::print(ostream & os, unsigned upper_precedence) const
 182 {
 183     debugmsg("idx print",LOGLEVEL_PRINT);
 184
 185     if (covariant) {
 186         os << "_";
 187     } else {
 188         os << "~";
 189     }
 190     if (symbolic) {
 191         os << name;
 192     } else {
 193         os << value;
 194     }
 195 }
 196
 197 bool idx::info(unsigned inf) const
 198 {
 199     if (inf==info_flags::idx) return true;
 200     return basic::info(inf);
 201 }
 202
 203 ex idx::subs(const lst & ls, const lst & lr) const
 204 {
 205     GINAC_ASSERT(ls.nops()==lr.nops());
 206 #ifdef DO_GINAC_ASSERT
 207     for (unsigned i=0; i<ls.nops(); i++) {
 208         GINAC_ASSERT(is_ex_exactly_of_type(ls.op(i),symbol)||
 209                is_ex_of_type(ls.op(i),idx));
 210     }
 211 #endif // def DO_GINAC_ASSERT
 212
 213     for (unsigned i=0; i<ls.nops(); i++) {
 214         if (is_equal(*(ls.op(i)).bp)) {
 215             return lr.op(i);
 216         }
 217     }
 218     return *this;
 219 }
 220
 221 // protected
 222
 223 int idx::compare_same_type(const basic & other) const
 224 {
 225     GINAC_ASSERT(is_of_type(other,idx));
 226     const idx & o=static_cast<const idx &>
 227                              (const_cast<basic &>(other));
 228
 229     if (covariant!=o.covariant) {
 230         // different co/contravariant
 231         return covariant ? -1 : 1;
 232     }
 233     if ((!symbolic) && (!o.symbolic)) {
 234         // non-symbolic, of equal type: compare values
 235         if (value==o.value) {
 236             return 0;
 237         }
 238         return value<o.value ? -1 : 1;
 239     }
 240     if (symbolic && o.symbolic) {
 241         // both symbolic: compare serials
 242         if (serial==o.serial) {
 243             return 0;
 244         }
 245         return serial<o.serial ? -1 : 1;
 246     }
 247     // one symbolic, one value: value is sorted first
 248     return o.symbolic ? -1 : 1;
 249 }
 250
 251 bool idx::is_equal_same_type(const basic & other) const
 252 {
 253     GINAC_ASSERT(is_of_type(other,idx));
 254     const idx & o=static_cast<const idx &>
 255                              (const_cast<basic &>(other));
 256
 257     if (covariant!=o.covariant) return false;
 258     if (symbolic!=o.symbolic) return false;
 259     if (symbolic && o.symbolic) return serial==o.serial;
 260     return value==o.value;
 261 }
 262
 263 unsigned idx::calchash(void) const
 264 {
 265     hashvalue=golden_ratio_hash(golden_ratio_hash(tinfo_key ^ serial));
 266     setflag(status_flags::hash_calculated);
 267     return hashvalue;
 268 }
 269
 270 //////////
 271 // new virtual functions which can be overridden by derived classes
 272 //////////
 273
 274 // public
 275
 276 bool idx::is_co_contra_pair(const basic & other) const
 277 {
 278     // like is_equal_same_type(), but tests for different covariant status
 279     GINAC_ASSERT(is_of_type(other,idx));
 280     const idx & o=static_cast<const idx &>
 281                              (const_cast<basic &>(other));
 282
 283     if (covariant==o.covariant) return false;
 284     if (symbolic!=o.symbolic) return false;
 285     if (symbolic && o.symbolic) return serial==o.serial;
 286     return value==o.value;
 287 }
 288
 289 bool idx::is_symbolic(void) const
 290 {
 291     return symbolic;
 292 }
 293
 294 unsigned idx::get_value(void) const
 295 {
 296     return value;
 297 }
 298
 299 bool idx::is_covariant(void) const
 300 {
 301     return covariant;
 302 }
 303
 304 ex idx::toggle_covariant(void) const
 305 {
 306     idx * i_copy=static_cast<idx *>(duplicate());
 307     i_copy->covariant = !i_copy->covariant;
 308     i_copy->clearflag(status_flags::hash_calculated);
 309     return i_copy->setflag(status_flags::dynallocated);
 310 }
 311
 312 //////////
 313 // non-virtual functions in this class
 314 //////////
 315
 316 // none
 317
 318 //////////
 319 // static member variables
 320 //////////
 321
 322 // protected
 323
 324 unsigned idx::next_serial=0;
 325
 326 //////////
 327 // global constants
 328 //////////
 329
 330 const idx some_idx;
 331 const type_info & typeid_idx=typeid(some_idx);
 332
 333 //////////
 334 // other functions
 335 //////////
 336
 337 int canonicalize_indices(exvector & iv, bool antisymmetric)
 338 {
 339     if (iv.size()<2) {
 340         // nothing do to for 0 or 1 indices
 341         return INT_MAX;
 342     }
 343
 344     bool something_changed=false;
 345     int sig=1;
 346     // simple bubble sort algorithm should be sufficient for the small number of indices needed
 347     exvector::const_iterator last_idx=iv.end();
 348     exvector::const_iterator next_to_last_idx=iv.end()-1;
 349     for (exvector::iterator it1=iv.begin(); it1!=next_to_last_idx; ++it1) {
 350         for (exvector::iterator it2=it1+1; it2!=last_idx; ++it2) {
 351             int cmpval=(*it1).compare(*it2);
 352             if (cmpval==1) {
 353                 iter_swap(it1,it2);
 354                 something_changed=true;
 355                 if (antisymmetric) sig=-sig;
 356             } else if ((cmpval==0) && antisymmetric) {
 357                 something_changed=true;
 358                 sig=0;
 359             }
 360         }
 361     }
 362     return something_changed ? sig : INT_MAX;
 363 }
 364
 365 exvector idx_intersect(const exvector & iv1, const exvector & iv2)
 366 {
 367     // build a vector of symbolic indices contained in iv1 and iv2 simultaneously
 368     // assumes (but does not test) that each index occurs at most twice
 369     exvector iv_intersect;
 370     for (exvector::const_iterator cit1=iv1.begin(); cit1!=iv1.end(); ++cit1) {
 371         GINAC_ASSERT(is_ex_of_type(*cit1,idx));
 372         if (ex_to_idx(*cit1).is_symbolic()) {
 373             for (exvector::const_iterator cit2=iv2.begin(); cit2!=iv2.end(); ++cit2) {
 374                 GINAC_ASSERT(is_ex_of_type(*cit2,idx));
 375                 if ((*cit1).is_equal(*cit2)) {
 376                     iv_intersect.push_back(*cit1);
 377                     break;
 378                 }
 379             }
 380         }
 381     }
 382     return iv_intersect;
 383 }
 384
 385 #define TEST_PERMUTATION(A,B,C,P) \
 386     if ((iv3[B].is_equal(iv2[0]))&&(iv3[C].is_equal(iv2[1]))) { \
 387         if (antisymmetric) *sig=P; \
 388         return iv3[A]; \
 389     }
 390
 391 ex permute_free_index_to_front(const exvector & iv3, const exvector & iv2,
 392                                bool antisymmetric, int * sig)
 393 {
 394     // match (return value,iv2) to iv3 by permuting indices
 395     // iv3 is always cyclic
 396
 397     GINAC_ASSERT(iv3.size()==3);
 398     GINAC_ASSERT(iv2.size()==2);
 399
 400     *sig=1;
 401
 402     TEST_PERMUTATION(0,1,2,  1);
 403     TEST_PERMUTATION(0,2,1, -1);
 404     TEST_PERMUTATION(1,0,2, -1);
 405     TEST_PERMUTATION(1,2,0,  1);
 406     TEST_PERMUTATION(2,0,1,  1);
 407     TEST_PERMUTATION(2,1,0, -1);
 408     throw(std::logic_error("permute_free_index_to_front(): no valid permutation found"));
 409 }
 410
 411 unsigned subs_index_in_exvector(exvector & v, const ex & is, const ex & ir)
 412 {
 413     exvector::iterator it;
 414     unsigned replacements=0;
 415     unsigned current_replacements;
 416
 417     GINAC_ASSERT(is_ex_of_type(is,idx));
 418     GINAC_ASSERT(is_ex_of_type(ir,idx));
 419
 420     for (it=v.begin(); it!=v.end(); ++it) {
 421         current_replacements=count_index(*it,is);
 422         if (current_replacements>0) {
 423             (*it)=(*it).subs(is==ir);
 424         }
 425         replacements += current_replacements;
 426     }
 427     return replacements;
 428 }
 429
 430 unsigned count_index(const ex & e, const ex & i)
 431 {
 432     exvector idxv=e.get_indices();
 433     unsigned count=0;
 434     for (exvector::const_iterator cit=idxv.begin(); cit!=idxv.end(); ++cit) {
 435         if ((*cit).is_equal(i)) count++;
 436     }
 437     return count;
 438 }
 439
 440 ex subs_indices(const ex & e, const exvector & idxv_subs,
 441                 const exvector & idxv_repl)
 442 {
 443     GINAC_ASSERT(idxv_subs.size()==idxv_repl.size());
 444     ex res=e;
 445     for (unsigned i=0; i<idxv_subs.size(); ++i) {
 446         res=res.subs(idxv_subs[i]==idxv_repl[i]);
 447     }
 448     return res;
 449 }
 450
 451 #ifndef NO_GINAC_NAMESPACE
 452 } // namespace GiNaC
 453 #endif // ndef NO_GINAC_NAMESPACE