ginac/idx.cpp

   1 /** @file idx.cpp
   2  *
   3  *  Implementation of GiNaC's indices. */
   4
   5 /*
   6  *  GiNaC Copyright (C) 1999-2001 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 <iostream>
  24 #include <stdexcept>
  25
  26 #include "idx.h"
  27 #include "symbol.h"
  28 #include "lst.h"
  29 #include "print.h"
  30 #include "archive.h"
  31 #include "utils.h"
  32
  33 namespace GiNaC {
  34
  35 GINAC_IMPLEMENT_REGISTERED_CLASS(idx, basic)
  36 GINAC_IMPLEMENT_REGISTERED_CLASS(varidx, idx)
  37 GINAC_IMPLEMENT_REGISTERED_CLASS(spinidx, varidx)
  38
  39 //////////
  40 // default ctor, dtor, copy ctor, assignment operator and helpers
  41 //////////
  42
  43 idx::idx() : inherited(TINFO_idx) {}
  44
  45 varidx::varidx() : covariant(false)
  46 {
  47         tinfo_key = TINFO_varidx;
  48 }
  49
  50 spinidx::spinidx() : dotted(false)
  51 {
  52         tinfo_key = TINFO_spinidx;
  53 }
  54
  55 void idx::copy(const idx & other)
  56 {
  57         inherited::copy(other);
  58         value = other.value;
  59         dim = other.dim;
  60 }
  61
  62 void varidx::copy(const varidx & other)
  63 {
  64         inherited::copy(other);
  65         covariant = other.covariant;
  66 }
  67
  68 void spinidx::copy(const spinidx & other)
  69 {
  70         inherited::copy(other);
  71         dotted = other.dotted;
  72 }
  73
  74 DEFAULT_DESTROY(idx)
  75 DEFAULT_DESTROY(varidx)
  76 DEFAULT_DESTROY(spinidx)
  77
  78 //////////
  79 // other constructors
  80 //////////
  81
  82 idx::idx(const ex & v, const ex & d) : inherited(TINFO_idx), value(v), dim(d)
  83 {
  84         if (is_dim_numeric())
  85                 if (!dim.info(info_flags::posint))
  86                         throw(std::invalid_argument("dimension of space must be a positive integer"));
  87 }
  88
  89 varidx::varidx(const ex & v, const ex & d, bool cov) : inherited(v, d), covariant(cov)
  90 {
  91         tinfo_key = TINFO_varidx;
  92 }
  93
  94 spinidx::spinidx(const ex & v, const ex & d, bool cov, bool dot) : inherited(v, d, cov), dotted(dot)
  95 {
  96         tinfo_key = TINFO_spinidx;
  97 }
  98
  99 //////////
 100 // archiving
 101 //////////
 102
 103 idx::idx(const archive_node &n, const lst &sym_lst) : inherited(n, sym_lst)
 104 {
 105         n.find_ex("value", value, sym_lst);
 106         n.find_ex("dim", dim, sym_lst);
 107 }
 108
 109 varidx::varidx(const archive_node &n, const lst &sym_lst) : inherited(n, sym_lst)
 110 {
 111         n.find_bool("covariant", covariant);
 112 }
 113
 114 spinidx::spinidx(const archive_node &n, const lst &sym_lst) : inherited(n, sym_lst)
 115 {
 116         n.find_bool("dotted", dotted);
 117 }
 118
 119 void idx::archive(archive_node &n) const
 120 {
 121         inherited::archive(n);
 122         n.add_ex("value", value);
 123         n.add_ex("dim", dim);
 124 }
 125
 126 void varidx::archive(archive_node &n) const
 127 {
 128         inherited::archive(n);
 129         n.add_bool("covariant", covariant);
 130 }
 131
 132 void spinidx::archive(archive_node &n) const
 133 {
 134         inherited::archive(n);
 135         n.add_bool("dotted", dotted);
 136 }
 137
 138 DEFAULT_UNARCHIVE(idx)
 139 DEFAULT_UNARCHIVE(varidx)
 140 DEFAULT_UNARCHIVE(spinidx)
 141
 142 //////////
 143 // functions overriding virtual functions from base classes
 144 //////////
 145
 146 void idx::print(const print_context & c, unsigned level) const
 147 {
 148         if (is_of_type(c, print_tree)) {
 149
 150                 c.s << std::string(level, ' ') << class_name()
 151                     << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
 152                     << std::endl;
 153                 unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
 154                 value.print(c, level + delta_indent);
 155                 dim.print(c, level + delta_indent);
 156
 157         } else {
 158
 159                 if (is_a<print_latex>(c))
 160                         c.s << "_{";
 161                 else
 162                         c.s << ".";
 163                 bool need_parens = !(is_ex_exactly_of_type(value, numeric) || is_ex_of_type(value, symbol));
 164                 if (need_parens)
 165                         c.s << "(";
 166                 value.print(c);
 167                 if (need_parens)
 168                         c.s << ")";
 169                 if (is_a<print_latex>(c))
 170                         c.s << "}";
 171         }
 172 }
 173
 174 void varidx::print(const print_context & c, unsigned level) const
 175 {
 176         if (is_of_type(c, print_tree)) {
 177
 178                 c.s << std::string(level, ' ') << class_name()
 179                     << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
 180                     << (covariant ? ", covariant" : ", contravariant")
 181                     << std::endl;
 182                 unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
 183                 value.print(c, level + delta_indent);
 184                 dim.print(c, level + delta_indent);
 185
 186         } else {
 187                 if (is_a<print_latex>(c)) {
 188                         if (covariant)
 189                                 c.s << "_{";
 190                         else
 191                                 c.s << "^{";
 192                 } else {
 193                         if (covariant)
 194                                 c.s << ".";
 195                         else
 196                                 c.s << "~";
 197                 }
 198                 bool need_parens = !(is_ex_exactly_of_type(value, numeric) || is_ex_of_type(value, symbol));
 199                 if (need_parens)
 200                         c.s << "(";
 201                 value.print(c);
 202                 if (need_parens)
 203                         c.s << ")";
 204                 if (is_a<print_latex>(c))
 205                         c.s << "}";
 206         }
 207 }
 208
 209 void spinidx::print(const print_context & c, unsigned level) const
 210 {
 211         if (is_of_type(c, print_tree)) {
 212
 213                 c.s << std::string(level, ' ') << class_name()
 214                     << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
 215                     << (covariant ? ", covariant" : ", contravariant")
 216                     << (dotted ? ", dotted" : ", undotted")
 217                     << std::endl;
 218                 unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
 219                 value.print(c, level + delta_indent);
 220                 dim.print(c, level + delta_indent);
 221
 222         } else {
 223
 224                 bool is_tex = is_of_type(c, print_latex);
 225                 if (is_tex) {
 226                         if (covariant)
 227                                 c.s << "_{";
 228                         else
 229                                 c.s << "^{";
 230                 } else {
 231                         if (covariant)
 232                                 c.s << ".";
 233                         else
 234                                 c.s << "~";
 235                 }
 236                 if (dotted) {
 237                         if (is_tex)
 238                                 c.s << "\\dot{";
 239                         else
 240                                 c.s << "*";
 241                 }
 242                 bool need_parens = !(is_ex_exactly_of_type(value, numeric) || is_ex_of_type(value, symbol));
 243                 if (need_parens)
 244                         c.s << "(";
 245                 value.print(c);
 246                 if (need_parens)
 247                         c.s << ")";
 248                 if (is_tex && dotted)
 249                         c.s << "}";
 250                 if (is_tex)
 251                         c.s << "}";
 252         }
 253 }
 254
 255 bool idx::info(unsigned inf) const
 256 {
 257         if (inf == info_flags::idx)
 258                 return true;
 259         return inherited::info(inf);
 260 }
 261
 262 unsigned idx::nops() const
 263 {
 264         // don't count the dimension as that is not really a sub-expression
 265         return 1;
 266 }
 267
 268 ex & idx::let_op(int i)
 269 {
 270         GINAC_ASSERT(i == 0);
 271         return value;
 272 }
 273
 274 /** Returns order relation between two indices of the same type. The order
 275  *  must be such that dummy indices lie next to each other. */
 276 int idx::compare_same_type(const basic & other) const
 277 {
 278         GINAC_ASSERT(is_a<idx>(other));
 279         const idx &o = static_cast<const idx &>(other);
 280
 281         int cmpval = value.compare(o.value);
 282         if (cmpval)
 283                 return cmpval;
 284         return dim.compare(o.dim);
 285 }
 286
 287 bool idx::match_same_type(const basic & other) const
 288 {
 289         GINAC_ASSERT(is_a<idx>(other));
 290         const idx &o = static_cast<const idx &>(other);
 291
 292         return dim.is_equal(o.dim);
 293 }
 294
 295 int varidx::compare_same_type(const basic & other) const
 296 {
 297         GINAC_ASSERT(is_a<varidx>(other));
 298         const varidx &o = static_cast<const varidx &>(other);
 299
 300         int cmpval = inherited::compare_same_type(other);
 301         if (cmpval)
 302                 return cmpval;
 303
 304         // Check variance last so dummy indices will end up next to each other
 305         if (covariant != o.covariant)
 306                 return covariant ? -1 : 1;
 307         return 0;
 308 }
 309
 310 bool varidx::match_same_type(const basic & other) const
 311 {
 312         GINAC_ASSERT(is_a<varidx>(other));
 313         const varidx &o = static_cast<const varidx &>(other);
 314
 315         if (covariant != o.covariant)
 316                 return false;
 317         return inherited::match_same_type(other);
 318 }
 319
 320 int spinidx::compare_same_type(const basic & other) const
 321 {
 322         GINAC_ASSERT(is_a<spinidx>(other));
 323         const spinidx &o = static_cast<const spinidx &>(other);
 324
 325         // Check dottedness first so dummy indices will end up next to each other
 326         if (dotted != o.dotted)
 327                 return dotted ? -1 : 1;
 328
 329         int cmpval = inherited::compare_same_type(other);
 330         if (cmpval)
 331                 return cmpval;
 332
 333         return 0;
 334 }
 335
 336 bool spinidx::match_same_type(const basic & other) const
 337 {
 338         GINAC_ASSERT(is_a<spinidx>(other));
 339         const spinidx &o = static_cast<const spinidx &>(other);
 340
 341         if (dotted != o.dotted)
 342                 return false;
 343         return inherited::match_same_type(other);
 344 }
 345
 346 /** By default, basic::evalf would evaluate the index value but we don't want
 347  *  a.1 to become a.(1.0). */
 348 ex idx::evalf(int level) const
 349 {
 350         return *this;
 351 }
 352
 353 ex idx::subs(const lst & ls, const lst & lr, bool no_pattern) const
 354 {
 355         GINAC_ASSERT(ls.nops() == lr.nops());
 356
 357         // First look for index substitutions
 358         for (unsigned i=0; i<ls.nops(); i++) {
 359                 if (is_equal(ex_to<basic>(ls.op(i)))) {
 360
 361                         // Substitution index->index
 362                         if (is_ex_of_type(lr.op(i), idx))
 363                                 return lr.op(i);
 364
 365                         // Otherwise substitute value
 366                         idx *i_copy = static_cast<idx *>(duplicate());
 367                         i_copy->value = lr.op(i);
 368                         i_copy->clearflag(status_flags::hash_calculated);
 369                         return i_copy->setflag(status_flags::dynallocated);
 370                 }
 371         }
 372
 373         // None, substitute objects in value (not in dimension)
 374         const ex &subsed_value = value.subs(ls, lr, no_pattern);
 375         if (are_ex_trivially_equal(value, subsed_value))
 376                 return *this;
 377
 378         idx *i_copy = static_cast<idx *>(duplicate());
 379         i_copy->value = subsed_value;
 380         i_copy->clearflag(status_flags::hash_calculated);
 381         return i_copy->setflag(status_flags::dynallocated);
 382 }
 383
 384 /** Implementation of ex::diff() for an index always returns 0.
 385  *
 386  *  @see ex::diff */
 387 ex idx::derivative(const symbol & s) const
 388 {
 389         return _ex0;
 390 }
 391
 392 //////////
 393 // new virtual functions
 394 //////////
 395
 396 bool idx::is_dummy_pair_same_type(const basic & other) const
 397 {
 398         const idx &o = static_cast<const idx &>(other);
 399
 400         // Only pure symbols form dummy pairs, "2n+1" doesn't
 401         if (!is_ex_of_type(value, symbol))
 402                 return false;
 403
 404         // Value must be equal, of course
 405         if (!value.is_equal(o.value))
 406                 return false;
 407
 408         // Also the dimension
 409         return dim.is_equal(o.dim);
 410 }
 411
 412 bool varidx::is_dummy_pair_same_type(const basic & other) const
 413 {
 414         const varidx &o = static_cast<const varidx &>(other);
 415
 416         // Variance must be opposite
 417         if (covariant == o.covariant)
 418                 return false;
 419
 420         return inherited::is_dummy_pair_same_type(other);
 421 }
 422
 423 bool spinidx::is_dummy_pair_same_type(const basic & other) const
 424 {
 425         const spinidx &o = static_cast<const spinidx &>(other);
 426
 427         // Dottedness must be the same
 428         if (dotted != o.dotted)
 429                 return false;
 430
 431         return inherited::is_dummy_pair_same_type(other);
 432 }
 433
 434
 435 //////////
 436 // non-virtual functions
 437 //////////
 438
 439 ex varidx::toggle_variance(void) const
 440 {
 441         varidx *i_copy = static_cast<varidx *>(duplicate());
 442         i_copy->covariant = !i_copy->covariant;
 443         i_copy->clearflag(status_flags::hash_calculated);
 444         return i_copy->setflag(status_flags::dynallocated);
 445 }
 446
 447 ex spinidx::toggle_dot(void) const
 448 {
 449         spinidx *i_copy = static_cast<spinidx *>(duplicate());
 450         i_copy->dotted = !i_copy->dotted;
 451         i_copy->clearflag(status_flags::hash_calculated);
 452         return i_copy->setflag(status_flags::dynallocated);
 453 }
 454
 455 ex spinidx::toggle_variance_dot(void) const
 456 {
 457         spinidx *i_copy = static_cast<spinidx *>(duplicate());
 458         i_copy->covariant = !i_copy->covariant;
 459         i_copy->dotted = !i_copy->dotted;
 460         i_copy->clearflag(status_flags::hash_calculated);
 461         return i_copy->setflag(status_flags::dynallocated);
 462 }
 463
 464 //////////
 465 // global functions
 466 //////////
 467
 468 bool is_dummy_pair(const idx & i1, const idx & i2)
 469 {
 470         // The indices must be of exactly the same type
 471         if (i1.tinfo() != i2.tinfo())
 472                 return false;
 473
 474         // Same type, let the indices decide whether they are paired
 475         return i1.is_dummy_pair_same_type(i2);
 476 }
 477
 478 bool is_dummy_pair(const ex & e1, const ex & e2)
 479 {
 480         // The expressions must be indices
 481         if (!is_ex_of_type(e1, idx) || !is_ex_of_type(e2, idx))
 482                 return false;
 483
 484         return is_dummy_pair(ex_to<idx>(e1), ex_to<idx>(e2));
 485 }
 486
 487 void find_free_and_dummy(exvector::const_iterator it, exvector::const_iterator itend, exvector & out_free, exvector & out_dummy)
 488 {
 489         out_free.clear();
 490         out_dummy.clear();
 491
 492         // No indices? Then do nothing
 493         if (it == itend)
 494                 return;
 495
 496         // Only one index? Then it is a free one if it's not numeric
 497         if (itend - it == 1) {
 498                 if (ex_to<idx>(*it).is_symbolic())
 499                         out_free.push_back(*it);
 500                 return;
 501         }
 502
 503         // Sort index vector. This will cause dummy indices come to lie next
 504         // to each other (because the sort order is defined to guarantee this).
 505         exvector v(it, itend);
 506         shaker_sort(v.begin(), v.end(), ex_is_less(), ex_swap());
 507
 508         // Find dummy pairs and free indices
 509         it = v.begin(); itend = v.end();
 510         exvector::const_iterator last = it++;
 511         while (it != itend) {
 512                 if (is_dummy_pair(*it, *last)) {
 513                         out_dummy.push_back(*last);
 514                         it++;
 515                         if (it == itend)
 516                                 return;
 517                 } else {
 518                         if (!it->is_equal(*last) && ex_to<idx>(*last).is_symbolic())
 519                                 out_free.push_back(*last);
 520                 }
 521                 last = it++;
 522         }
 523         if (ex_to<idx>(*last).is_symbolic())
 524                 out_free.push_back(*last);
 525 }
 526
 527 } // namespace GiNaC