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