ginac/idx.cpp

   1 /** @file idx.cpp
   2  *
   3  *  Implementation of GiNaC's indices. */
   4
   5 /*
   6  *  GiNaC Copyright (C) 1999-2002 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                         c.s << "{";
 189                 else {
 190                         if (covariant)
 191                                 c.s << ".";
 192                         else
 193                                 c.s << "~";
 194                 }
 195                 bool need_parens = !(is_ex_exactly_of_type(value, numeric) || is_ex_of_type(value, symbol));
 196                 if (need_parens)
 197                         c.s << "(";
 198                 value.print(c);
 199                 if (need_parens)
 200                         c.s << ")";
 201                 if (is_a<print_latex>(c))
 202                         c.s << "}";
 203         }
 204 }
 205
 206 void spinidx::print(const print_context & c, unsigned level) const
 207 {
 208         if (is_of_type(c, print_tree)) {
 209
 210                 c.s << std::string(level, ' ') << class_name()
 211                     << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
 212                     << (covariant ? ", covariant" : ", contravariant")
 213                     << (dotted ? ", dotted" : ", undotted")
 214                     << std::endl;
 215                 unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
 216                 value.print(c, level + delta_indent);
 217                 dim.print(c, level + delta_indent);
 218
 219         } else {
 220
 221                 bool is_tex = is_of_type(c, print_latex);
 222                 if (is_tex) {
 223                         if (covariant)
 224                                 c.s << "_{";
 225                         else
 226                                 c.s << "^{";
 227                 } else {
 228                         if (covariant)
 229                                 c.s << ".";
 230                         else
 231                                 c.s << "~";
 232                 }
 233                 if (dotted) {
 234                         if (is_tex)
 235                                 c.s << "\\dot{";
 236                         else
 237                                 c.s << "*";
 238                 }
 239                 bool need_parens = !(is_ex_exactly_of_type(value, numeric) || is_ex_of_type(value, symbol));
 240                 if (need_parens)
 241                         c.s << "(";
 242                 value.print(c);
 243                 if (need_parens)
 244                         c.s << ")";
 245                 if (is_tex && dotted)
 246                         c.s << "}";
 247                 if (is_tex)
 248                         c.s << "}";
 249         }
 250 }
 251
 252 bool idx::info(unsigned inf) const
 253 {
 254         if (inf == info_flags::idx)
 255                 return true;
 256         return inherited::info(inf);
 257 }
 258
 259 unsigned idx::nops() const
 260 {
 261         // don't count the dimension as that is not really a sub-expression
 262         return 1;
 263 }
 264
 265 ex & idx::let_op(int i)
 266 {
 267         GINAC_ASSERT(i == 0);
 268         return value;
 269 }
 270
 271 /** Returns order relation between two indices of the same type. The order
 272  *  must be such that dummy indices lie next to each other. */
 273 int idx::compare_same_type(const basic & other) const
 274 {
 275         GINAC_ASSERT(is_a<idx>(other));
 276         const idx &o = static_cast<const idx &>(other);
 277
 278         int cmpval = value.compare(o.value);
 279         if (cmpval)
 280                 return cmpval;
 281         return dim.compare(o.dim);
 282 }
 283
 284 bool idx::match_same_type(const basic & other) const
 285 {
 286         GINAC_ASSERT(is_a<idx>(other));
 287         const idx &o = static_cast<const idx &>(other);
 288
 289         return dim.is_equal(o.dim);
 290 }
 291
 292 int varidx::compare_same_type(const basic & other) const
 293 {
 294         GINAC_ASSERT(is_a<varidx>(other));
 295         const varidx &o = static_cast<const varidx &>(other);
 296
 297         int cmpval = inherited::compare_same_type(other);
 298         if (cmpval)
 299                 return cmpval;
 300
 301         // Check variance last so dummy indices will end up next to each other
 302         if (covariant != o.covariant)
 303                 return covariant ? -1 : 1;
 304         return 0;
 305 }
 306
 307 bool varidx::match_same_type(const basic & other) const
 308 {
 309         GINAC_ASSERT(is_a<varidx>(other));
 310         const varidx &o = static_cast<const varidx &>(other);
 311
 312         if (covariant != o.covariant)
 313                 return false;
 314         return inherited::match_same_type(other);
 315 }
 316
 317 int spinidx::compare_same_type(const basic & other) const
 318 {
 319         GINAC_ASSERT(is_a<spinidx>(other));
 320         const spinidx &o = static_cast<const spinidx &>(other);
 321
 322         // Check dottedness first so dummy indices will end up next to each other
 323         if (dotted != o.dotted)
 324                 return dotted ? -1 : 1;
 325
 326         int cmpval = inherited::compare_same_type(other);
 327         if (cmpval)
 328                 return cmpval;
 329
 330         return 0;
 331 }
 332
 333 bool spinidx::match_same_type(const basic & other) const
 334 {
 335         GINAC_ASSERT(is_a<spinidx>(other));
 336         const spinidx &o = static_cast<const spinidx &>(other);
 337
 338         if (dotted != o.dotted)
 339                 return false;
 340         return inherited::match_same_type(other);
 341 }
 342
 343 /** By default, basic::evalf would evaluate the index value but we don't want
 344  *  a.1 to become a.(1.0). */
 345 ex idx::evalf(int level) const
 346 {
 347         return *this;
 348 }
 349
 350 ex idx::subs(const lst & ls, const lst & lr, bool no_pattern) const
 351 {
 352         GINAC_ASSERT(ls.nops() == lr.nops());
 353
 354         // First look for index substitutions
 355         for (unsigned i=0; i<ls.nops(); i++) {
 356                 if (is_equal(ex_to<basic>(ls.op(i)))) {
 357
 358                         // Substitution index->index
 359                         if (is_ex_of_type(lr.op(i), idx))
 360                                 return lr.op(i);
 361
 362                         // Otherwise substitute value
 363                         idx *i_copy = static_cast<idx *>(duplicate());
 364                         i_copy->value = lr.op(i);
 365                         i_copy->clearflag(status_flags::hash_calculated);
 366                         return i_copy->setflag(status_flags::dynallocated);
 367                 }
 368         }
 369
 370         // None, substitute objects in value (not in dimension)
 371         const ex &subsed_value = value.subs(ls, lr, no_pattern);
 372         if (are_ex_trivially_equal(value, subsed_value))
 373                 return *this;
 374
 375         idx *i_copy = static_cast<idx *>(duplicate());
 376         i_copy->value = subsed_value;
 377         i_copy->clearflag(status_flags::hash_calculated);
 378         return i_copy->setflag(status_flags::dynallocated);
 379 }
 380
 381 /** Implementation of ex::diff() for an index always returns 0.
 382  *
 383  *  @see ex::diff */
 384 ex idx::derivative(const symbol & s) const
 385 {
 386         return _ex0;
 387 }
 388
 389 //////////
 390 // new virtual functions
 391 //////////
 392
 393 bool idx::is_dummy_pair_same_type(const basic & other) const
 394 {
 395         const idx &o = static_cast<const idx &>(other);
 396
 397         // Only pure symbols form dummy pairs, "2n+1" doesn't
 398         if (!is_ex_of_type(value, symbol))
 399                 return false;
 400
 401         // Value must be equal, of course
 402         if (!value.is_equal(o.value))
 403                 return false;
 404
 405         // Also the dimension
 406         return dim.is_equal(o.dim);
 407 }
 408
 409 bool varidx::is_dummy_pair_same_type(const basic & other) const
 410 {
 411         const varidx &o = static_cast<const varidx &>(other);
 412
 413         // Variance must be opposite
 414         if (covariant == o.covariant)
 415                 return false;
 416
 417         return inherited::is_dummy_pair_same_type(other);
 418 }
 419
 420 bool spinidx::is_dummy_pair_same_type(const basic & other) const
 421 {
 422         const spinidx &o = static_cast<const spinidx &>(other);
 423
 424         // Dottedness must be the same
 425         if (dotted != o.dotted)
 426                 return false;
 427
 428         return inherited::is_dummy_pair_same_type(other);
 429 }
 430
 431
 432 //////////
 433 // non-virtual functions
 434 //////////
 435
 436 ex varidx::toggle_variance(void) const
 437 {
 438         varidx *i_copy = static_cast<varidx *>(duplicate());
 439         i_copy->covariant = !i_copy->covariant;
 440         i_copy->clearflag(status_flags::hash_calculated);
 441         return i_copy->setflag(status_flags::dynallocated);
 442 }
 443
 444 ex spinidx::toggle_dot(void) const
 445 {
 446         spinidx *i_copy = static_cast<spinidx *>(duplicate());
 447         i_copy->dotted = !i_copy->dotted;
 448         i_copy->clearflag(status_flags::hash_calculated);
 449         return i_copy->setflag(status_flags::dynallocated);
 450 }
 451
 452 ex spinidx::toggle_variance_dot(void) const
 453 {
 454         spinidx *i_copy = static_cast<spinidx *>(duplicate());
 455         i_copy->covariant = !i_copy->covariant;
 456         i_copy->dotted = !i_copy->dotted;
 457         i_copy->clearflag(status_flags::hash_calculated);
 458         return i_copy->setflag(status_flags::dynallocated);
 459 }
 460
 461 //////////
 462 // global functions
 463 //////////
 464
 465 bool is_dummy_pair(const idx & i1, const idx & i2)
 466 {
 467         // The indices must be of exactly the same type
 468         if (i1.tinfo() != i2.tinfo())
 469                 return false;
 470
 471         // Same type, let the indices decide whether they are paired
 472         return i1.is_dummy_pair_same_type(i2);
 473 }
 474
 475 bool is_dummy_pair(const ex & e1, const ex & e2)
 476 {
 477         // The expressions must be indices
 478         if (!is_ex_of_type(e1, idx) || !is_ex_of_type(e2, idx))
 479                 return false;
 480
 481         return is_dummy_pair(ex_to<idx>(e1), ex_to<idx>(e2));
 482 }
 483
 484 void find_free_and_dummy(exvector::const_iterator it, exvector::const_iterator itend, exvector & out_free, exvector & out_dummy)
 485 {
 486         out_free.clear();
 487         out_dummy.clear();
 488
 489         // No indices? Then do nothing
 490         if (it == itend)
 491                 return;
 492
 493         // Only one index? Then it is a free one if it's not numeric
 494         if (itend - it == 1) {
 495                 if (ex_to<idx>(*it).is_symbolic())
 496                         out_free.push_back(*it);
 497                 return;
 498         }
 499
 500         // Sort index vector. This will cause dummy indices come to lie next
 501         // to each other (because the sort order is defined to guarantee this).
 502         exvector v(it, itend);
 503         shaker_sort(v.begin(), v.end(), ex_is_less(), ex_swap());
 504
 505         // Find dummy pairs and free indices
 506         it = v.begin(); itend = v.end();
 507         exvector::const_iterator last = it++;
 508         while (it != itend) {
 509                 if (is_dummy_pair(*it, *last)) {
 510                         out_dummy.push_back(*last);
 511                         it++;
 512                         if (it == itend)
 513                                 return;
 514                 } else {
 515                         if (!it->is_equal(*last) && ex_to<idx>(*last).is_symbolic())
 516                                 out_free.push_back(*last);
 517                 }
 518                 last = it++;
 519         }
 520         if (ex_to<idx>(*last).is_symbolic())
 521                 out_free.push_back(*last);
 522 }
 523
 524 } // namespace GiNaC