]> www.ginac.de Git - ginac.git/blob - ginac/idx.cpp
9f31a9ad8f997788f484d274ca1a07059745e76a
[ginac.git] / 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 <stdexcept>
24
25 #include "idx.h"
26 #include "symbol.h"
27 #include "lst.h"
28 #include "print.h"
29 #include "archive.h"
30 #include "utils.h"
31 #include "debugmsg.h"
32
33 namespace GiNaC {
34
35 GINAC_IMPLEMENT_REGISTERED_CLASS(idx, basic)
36 GINAC_IMPLEMENT_REGISTERED_CLASS(varidx, idx)
37
38 //////////
39 // default constructor, destructor, copy constructor assignment operator and helpers
40 //////////
41
42 idx::idx() : inherited(TINFO_idx)
43 {
44         debugmsg("idx default constructor", LOGLEVEL_CONSTRUCT);
45 }
46
47 varidx::varidx() : covariant(false)
48 {
49         debugmsg("varidx default constructor", LOGLEVEL_CONSTRUCT);
50         tinfo_key = TINFO_varidx;
51 }
52
53 void idx::copy(const idx & other)
54 {
55         inherited::copy(other);
56         value = other.value;
57         dim = other.dim;
58 }
59
60 void varidx::copy(const varidx & other)
61 {
62         inherited::copy(other);
63         covariant = other.covariant;
64 }
65
66 DEFAULT_DESTROY(idx)
67 DEFAULT_DESTROY(varidx)
68
69 //////////
70 // other constructors
71 //////////
72
73 idx::idx(const ex & v, const ex & d) : inherited(TINFO_idx), value(v), dim(d)
74 {
75         debugmsg("idx constructor from ex,ex", LOGLEVEL_CONSTRUCT);
76         if (is_dim_numeric())
77                 if (!dim.info(info_flags::posint))
78                         throw(std::invalid_argument("dimension of space must be a positive integer"));
79 }
80
81 varidx::varidx(const ex & v, const ex & d, bool cov) : inherited(v, d), covariant(cov)
82 {
83         debugmsg("varidx constructor from ex,ex,bool", LOGLEVEL_CONSTRUCT);
84         tinfo_key = TINFO_varidx;
85 }
86
87 //////////
88 // archiving
89 //////////
90
91 idx::idx(const archive_node &n, const lst &sym_lst) : inherited(n, sym_lst)
92 {
93         debugmsg("idx constructor from archive_node", LOGLEVEL_CONSTRUCT);
94         n.find_ex("value", value, sym_lst);
95         n.find_ex("dim", dim, sym_lst);
96 }
97
98 varidx::varidx(const archive_node &n, const lst &sym_lst) : inherited(n, sym_lst)
99 {
100         debugmsg("varidx constructor from archive_node", LOGLEVEL_CONSTRUCT);
101         n.find_bool("covariant", covariant);
102 }
103
104 void idx::archive(archive_node &n) const
105 {
106         inherited::archive(n);
107         n.add_ex("value", value);
108         n.add_ex("dim", dim);
109 }
110
111 void varidx::archive(archive_node &n) const
112 {
113         inherited::archive(n);
114         n.add_bool("covariant", covariant);
115 }
116
117 DEFAULT_UNARCHIVE(idx)
118 DEFAULT_UNARCHIVE(varidx)
119
120 //////////
121 // functions overriding virtual functions from bases classes
122 //////////
123
124 void idx::print(const print_context & c, unsigned level) const
125 {
126         debugmsg("idx print", LOGLEVEL_PRINT);
127
128         if (is_of_type(c, print_tree)) {
129
130                 c.s << std::string(level, ' ') << class_name()
131                     << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
132                     << std::endl;
133                 unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
134                 value.print(c, level + delta_indent);
135                 dim.print(c, level + delta_indent);
136
137         } else {
138
139                 c.s << ".";
140                 bool need_parens = !(is_ex_exactly_of_type(value, numeric) || is_ex_of_type(value, symbol));
141                 if (need_parens)
142                         c.s << "(";
143                 c.s << value;
144                 if (need_parens)
145                         c.s << ")";
146         }
147 }
148
149 void varidx::print(const print_context & c, unsigned level) const
150 {
151         debugmsg("varidx print", LOGLEVEL_PRINT);
152
153         if (is_of_type(c, print_tree)) {
154
155                 c.s << std::string(level, ' ') << class_name()
156                     << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
157                     << (covariant ? ", covariant" : ", contravariant")
158                     << std::endl;
159                 unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
160                 value.print(c, level + delta_indent);
161                 dim.print(c, level + delta_indent);
162
163         } else {
164
165                 if (covariant)
166                         c.s << ".";
167                 else
168                         c.s << "~";
169                 bool need_parens = !(is_ex_exactly_of_type(value, numeric) || is_ex_of_type(value, symbol));
170                 if (need_parens)
171                         c.s << "(";
172                 c.s << value;
173                 if (need_parens)
174                         c.s << ")";
175         }
176 }
177
178 bool idx::info(unsigned inf) const
179 {
180         if (inf == info_flags::idx)
181                 return true;
182         return inherited::info(inf);
183 }
184
185 unsigned idx::nops() const
186 {
187         // don't count the dimension as that is not really a sub-expression
188         return 1;
189 }
190
191 ex & idx::let_op(int i)
192 {
193         GINAC_ASSERT(i == 0);
194         return value;
195 }
196
197 /** Returns order relation between two indices of the same type. The order
198  *  must be such that dummy indices lie next to each other. */
199 int idx::compare_same_type(const basic & other) const
200 {
201         GINAC_ASSERT(is_of_type(other, idx));
202         const idx &o = static_cast<const idx &>(other);
203
204         int cmpval = value.compare(o.value);
205         if (cmpval)
206                 return cmpval;
207         return dim.compare(o.dim);
208 }
209
210 int varidx::compare_same_type(const basic & other) const
211 {
212         GINAC_ASSERT(is_of_type(other, varidx));
213         const varidx &o = static_cast<const varidx &>(other);
214
215         int cmpval = inherited::compare_same_type(other);
216         if (cmpval)
217                 return cmpval;
218
219         // Check variance last so dummy indices will end up next to each other
220         if (covariant != o.covariant)
221                 return covariant ? -1 : 1;
222         return 0;
223 }
224
225 ex idx::subs(const lst & ls, const lst & lr) const
226 {
227         GINAC_ASSERT(ls.nops() == lr.nops());
228
229         // First look for index substitutions
230         for (unsigned i=0; i<ls.nops(); i++) {
231                 if (is_equal(*(ls.op(i)).bp)) {
232
233                         // Substitution index->index
234                         if (is_ex_of_type(lr.op(i), idx))
235                                 return lr.op(i);
236
237                         // Otherwise substitute value
238                         idx *i_copy = static_cast<idx *>(duplicate());
239                         i_copy->value = lr.op(i);
240                         i_copy->clearflag(status_flags::hash_calculated);
241                         return i_copy->setflag(status_flags::dynallocated);
242                 }
243         }
244
245         // None, substitute objects in value (not in dimension)
246         const ex &subsed_value = value.subs(ls, lr);
247         if (are_ex_trivially_equal(value, subsed_value))
248                 return *this;
249
250         idx *i_copy = static_cast<idx *>(duplicate());
251         i_copy->value = subsed_value;
252         i_copy->clearflag(status_flags::hash_calculated);
253         return i_copy->setflag(status_flags::dynallocated);
254 }
255
256 //////////
257 // new virtual functions
258 //////////
259
260 bool idx::is_dummy_pair_same_type(const basic & other) const
261 {
262         const idx &o = static_cast<const idx &>(other);
263
264         // Only pure symbols form dummy pairs, "2n+1" doesn't
265         if (!is_ex_of_type(value, symbol))
266                 return false;
267
268         // Value must be equal, of course
269         if (!value.is_equal(o.value))
270                 return false;
271
272         // Also the dimension
273         return dim.is_equal(o.dim);
274 }
275
276 bool varidx::is_dummy_pair_same_type(const basic & other) const
277 {
278         const varidx &o = static_cast<const varidx &>(other);
279
280         // Variance must be opposite
281         if (covariant == o.covariant)
282                 return false;
283
284         return inherited::is_dummy_pair_same_type(other);
285 }
286
287 //////////
288 // non-virtual functions
289 //////////
290
291 ex varidx::toggle_variance(void) const
292 {
293         varidx *i_copy = static_cast<varidx *>(duplicate());
294         i_copy->covariant = !i_copy->covariant;
295         i_copy->clearflag(status_flags::hash_calculated);
296         return i_copy->setflag(status_flags::dynallocated);
297 }
298
299 //////////
300 // global functions
301 //////////
302
303 bool is_dummy_pair(const idx & i1, const idx & i2)
304 {
305         // The indices must be of exactly the same type
306         if (i1.tinfo() != i2.tinfo())
307                 return false;
308
309         // Same type, let the indices decide whether they are paired
310         return i1.is_dummy_pair_same_type(i2);
311 }
312
313 bool is_dummy_pair(const ex & e1, const ex & e2)
314 {
315         // The expressions must be indices
316         if (!is_ex_of_type(e1, idx) || !is_ex_of_type(e2, idx))
317                 return false;
318
319         return is_dummy_pair(ex_to_idx(e1), ex_to_idx(e2));
320 }
321
322 /** Bring a vector of indices into a canonic order. Dummy indices will lie
323  *  next to each other after the sorting. */
324 static void sort_index_vector(exvector &v)
325 {
326         // Nothing to sort if less than 2 elements
327         if (v.size() < 2)
328                 return;
329
330         // Simple bubble sort algorithm should be sufficient for the small
331         // number of indices expected
332         exvector::iterator it1 = v.begin(), itend = v.end(), next_to_last_idx = itend - 1;
333         while (it1 != next_to_last_idx) {
334                 exvector::iterator it2 = it1 + 1;
335                 while (it2 != itend) {
336                         if (it1->compare(*it2) > 0)
337                                 it1->swap(*it2);
338                         it2++;
339                 }
340                 it1++;
341         }
342 }
343
344
345 void find_free_and_dummy(exvector::const_iterator it, exvector::const_iterator itend, exvector & out_free, exvector & out_dummy)
346 {
347         out_free.clear();
348         out_dummy.clear();
349
350         // No indices? Then do nothing
351         if (it == itend)
352                 return;
353
354         // Only one index? Then it is a free one if it's not numeric
355         if (itend - it == 1) {
356                 if (ex_to_idx(*it).is_symbolic())
357                         out_free.push_back(*it);
358                 return;
359         }
360
361         // Sort index vector. This will cause dummy indices come to lie next
362         // to each other (because the sort order is defined to guarantee this).
363         exvector v(it, itend);
364         sort_index_vector(v);
365
366         // Find dummy pairs and free indices
367         it = v.begin(); itend = v.end();
368         exvector::const_iterator last = it++;
369         while (it != itend) {
370                 if (is_dummy_pair(*it, *last)) {
371                         out_dummy.push_back(*last);
372                         it++;
373                         if (it == itend)
374                                 return;
375                 } else {
376                         if (!it->is_equal(*last) && ex_to_idx(*last).is_symbolic())
377                                 out_free.push_back(*last);
378                 }
379                 last = it++;
380         }
381         if (ex_to_idx(*last).is_symbolic())
382                 out_free.push_back(*last);
383 }
384
385 exvector index_set_difference(const exvector & set1, const exvector & set2)
386 {
387         exvector ret;
388
389         exvector::const_iterator ait = set1.begin(), aitend = set1.end();
390         while (ait != aitend) {
391                 exvector::const_iterator bit = set2.begin(), bitend = set2.end();
392                 bool found = false;
393                 while (bit != bitend) {
394                         if (ait->is_equal(*bit)) {
395                                 found = true;
396                                 break;
397                         }
398                         bit++;
399                 }
400                 if (!found)
401                         ret.push_back(*ait);
402                 ait++;
403         }
404
405         return ret;
406 }
407
408 } // namespace GiNaC