glucat 0.12.0
index_set_imp.h
Go to the documentation of this file.
1#ifndef _GLUCAT_INDEX_SET_IMP_H
2#define _GLUCAT_INDEX_SET_IMP_H
3/***************************************************************************
4 GluCat : Generic library of universal Clifford algebra templates
5 index_set_imp.h : Implement a class for a set of non-zero integer indices
6 -------------------
7 begin : Sun 2001-12-09
8 copyright : (C) 2001-2016 by Paul C. Leopardi
9 ***************************************************************************
10
11 This library is free software: you can redistribute it and/or modify
12 it under the terms of the GNU Lesser General Public License as published
13 by the Free Software Foundation, either version 3 of the License, or
14 (at your option) any later version.
15
16 This library is distributed in the hope that it will be useful,
17 but WITHOUT ANY WARRANTY; without even the implied warranty of
18 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 GNU Lesser General Public License for more details.
20
21 You should have received a copy of the GNU Lesser General Public License
22 along with this library. If not, see <http://www.gnu.org/licenses/>.
23
24 ***************************************************************************
25 This library is based on a prototype written by Arvind Raja and was
26 licensed under the LGPL with permission of the author. See Arvind Raja,
27 "Object-oriented implementations of Clifford algebras in C++: a prototype",
28 in Ablamowicz, Lounesto and Parra (eds.)
29 "Clifford algebras with numeric and symbolic computations", Birkhauser, 1996.
30 ***************************************************************************
31 See also Arvind Raja's original header comments in glucat.h
32 ***************************************************************************/
33
34#include "glucat/index_set.h"
35
36#include <string>
37#include <sstream>
38
39namespace glucat
40{
41 // References for algorithms:
42 // [JA]: Joerg Arndt, "Algorithms for programmers", http://www.jjj.de/fxt/fxtbook.pdf
43 // Chapter 1, Bit wizardry, http://www.jjj.de/bitwizardry/bitwizardrypage.html
44 // [L]: Pertti Lounesto, "Clifford algebras and spinors", Cambridge UP, 1997.
45
46 template<const index_t LO, const index_t HI>
47 inline
48 auto
50 classname() -> const std::string
51 { return "index_set"; }
52
54 template<const index_t LO, const index_t HI>
56 index_set(const index_t idx)
57 { this->set(idx); }
58
60 template<const index_t LO, const index_t HI>
62 index_set(const bitset_t bst):
63 bitset_t(bst)
64 { }
65
67 template<const index_t LO, const index_t HI>
69 index_set(const set_value_t folded_val, const index_set_t frm, const bool prechecked)
70 {
71 if (!prechecked && folded_val >= (set_value_t(1) << frm.count()))
72 throw error_t("index_set(val,frm): cannot create: value gives an index set outside of frame");
73 const index_set_t folded_frame = frm.fold();
74 const index_t min_index = folded_frame.min();
75 const index_t skip = min_index > 0 ? 1 : 0;
76 const index_set_t folded_set = index_set_t(bitset_t(folded_val) << (min_index - skip - LO));
77 *this = folded_set.unfold(frm);
78 }
79
81 template<const index_t LO, const index_t HI>
83 index_set(const index_pair_t& range, const bool prechecked)
84 {
85 if (!prechecked && (range.first < LO || range.second > HI))
86 throw error_t("index_set(range): cannot create: range is too large");
87 const index_t begin_bit = (range.first < 0)
88 ? range.first-LO
89 : range.first-LO-1;
90 const index_t end_bit = (range.second < 0)
91 ? range.second-LO+1
92 : range.second-LO;
93 unsigned long mask = ( (end_bit == _GLUCAT_BITS_PER_ULONG)
94 ? -1UL
95 : (1UL << end_bit)-1UL)
96 & ~((1UL << begin_bit)-1UL);
97 *this = bitset_t(mask);
98 }
99
100 /// Constructor from string
101 template<const index_t LO, const index_t HI>
103 index_set(const std::string& str)
104 {
105 std::istringstream ss(str);
106 ss >> *this;
107 if (!ss)
108 throw error_t("index_set_t(str): could not parse string");
109 // Peek to see if the end of the string has been reached.
110 ss.peek();
111 if (!ss.eof())
112 throw error_t("index_set_t(str): could not parse entire string");
114
115 /// Equality
116 template<const index_t LO, const index_t HI>
117 inline
118 auto
120 operator== (const index_set_t rhs) const -> bool
122 const auto* pthis = static_cast<const bitset_t*>(this);
123 return *pthis == static_cast<bitset_t>(rhs);
124 }
127 template<const index_t LO, const index_t HI>
128 inline
129 auto
131 operator!= (const index_set_t rhs) const -> bool
132 {
133 const auto* pthis = static_cast<const bitset_t*>(this);
134 return *pthis != static_cast<bitset_t>(rhs);
136
137 /// Set complement: not
138 template<const index_t LO, const index_t HI>
139 inline
140 auto
142 operator~ () const -> index_set_t
143 { return bitset_t::operator~(); }
144
146 template<const index_t LO, const index_t HI>
147 inline
148 auto
152 bitset_t* pthis = this;
153 *pthis ^= static_cast<bitset_t>(rhs);
154 return *this;
156
157 /// Symmetric set difference: exclusive or
158 template<const index_t LO, const index_t HI>
159 inline
160 auto
162 const index_set<LO,HI>& rhs) -> const
166 using bitset_t = typename index_set_t::bitset_t;
167 return static_cast<bitset_t>(lhs) ^ static_cast<bitset_t>(rhs);
168 }
169
171 template<const index_t LO, const index_t HI>
172 inline
173 auto
176 {
177 bitset_t* pthis = this;
178 *pthis &= static_cast<bitset_t>(rhs);
179 return *this;
180 }
181
183 template<const index_t LO, const index_t HI>
184 inline
185 auto
187 const index_set<LO,HI>& rhs) -> const
189 {
191 using bitset_t = typename index_set_t::bitset_t;
192 return static_cast<bitset_t>(lhs) & static_cast<bitset_t>(rhs);
193 }
196 template<const index_t LO, const index_t HI>
197 inline
198 auto
201 {
202 bitset_t* pthis = this;
203 *pthis |= static_cast<bitset_t>(rhs);
204 return *this;
205 }
208 template<const index_t LO, const index_t HI>
209 inline
210 auto
212 const index_set<LO,HI>& rhs) -> const
214 {
215 using index_set_t = index_set<LO, HI>;
216 using bitset_t = typename index_set_t::bitset_t;
217 return static_cast<bitset_t>(lhs) | static_cast<bitset_t>(rhs);
218 }
219
221 template<const index_t LO, const index_t HI>
222 inline
223 auto
225 operator[] (const index_t idx) -> reference
226 { return reference(*this, idx); }
227
229 template<const index_t LO, const index_t HI>
230 inline
231 auto
233 operator[] (const index_t idx) const -> bool
234 { return this->test(idx); }
235
237 template<const index_t LO, const index_t HI>
238 inline
239 auto
241 test(const index_t idx) const -> bool
242 {
243 // Reference: [JA], 1.2.1
244 return (idx < 0)
245 ? bool(bitset_t::to_ulong() & (1UL << (idx - LO)))
246 : (idx > 0)
247 ? bool(bitset_t::to_ulong() & (1UL << (idx - LO - 1)))
248 : false;
249 }
250
252 template<const index_t LO, const index_t HI>
253 inline
254 auto
256 set() -> index_set_t&
257 {
258 bitset_t::set();
259 return *this;
260 }
261
263 template<const index_t LO, const index_t HI>
264 inline
265 auto
267 set(index_t idx) -> index_set_t&
268 {
269 if (idx > 0)
270 bitset_t::set(idx-LO-1);
271 else if (idx < 0)
272 bitset_t::set(idx-LO);
273 return *this;
274 }
275
277 template<const index_t LO, const index_t HI>
278 inline
279 auto
281 set(const index_t idx, const int val) -> index_set_t&
282 {
283 if (idx > 0)
284 bitset_t::set(idx-LO-1, val);
285 else if (idx < 0)
286 bitset_t::set(idx-LO, val);
287 return *this;
288 }
289
291 template<const index_t LO, const index_t HI>
292 inline
293 auto
296 {
297 bitset_t::reset();
298 return *this;
299 }
300
302 template<const index_t LO, const index_t HI>
303 inline
304 auto
306 reset(const index_t idx) -> index_set_t&
307 {
308 if (idx > 0)
309 bitset_t::reset(idx-LO-1);
310 else if (idx < 0)
311 bitset_t::reset(idx-LO);
312 return *this;
313 }
314
316 template<const index_t LO, const index_t HI>
317 inline
318 auto
321 {
322 bitset_t::flip();
323 return *this;
324 }
325
327 template<const index_t LO, const index_t HI>
328 inline
329 auto
331 flip(const index_t idx) -> index_set_t&
332 {
333 if (idx > 0)
334 bitset_t::flip(idx-LO-1);
335 else if (idx < 0)
336 bitset_t::flip(idx-LO);
337 return *this;
338 }
339
341 template<const index_t LO, const index_t HI>
342 inline
343 auto
345 count() const -> index_t
346 {
347 unsigned long val = bitset_t::to_ulong();
348 // Reference: [JA], 1.3
349 if (val == 0)
350 return 0;
351 else
352 {
353 index_t result = 1;
354 while (val &= val-1)
355 ++result;
356 return result;
357 }
358 }
359
361 template<const index_t LO, const index_t HI>
362 inline
363 auto
365 count_neg() const -> index_t
366 {
367 static const index_set_t lo_mask = bitset_t((1UL << -LO) - 1UL);
368 const index_set_t neg_part = *this & lo_mask;
369 return neg_part.count();
370 }
371
373 template<const index_t LO, const index_t HI>
374 inline
375 auto
377 count_pos() const -> index_t
378 {
379 const auto* pthis = static_cast<const bitset_t*>(this);
380 const index_set_t pos_part = *pthis >> -LO;
381 return pos_part.count();
382 }
383
384#if (_GLUCAT_BITS_PER_ULONG == 64)
386 template<const index_t LO, const index_t HI>
387 inline
388 auto
390 min() const -> index_t
391 {
392 // Reference: [JA], 1.3
393 unsigned long val = bitset_t::to_ulong();
394 if (val == 0)
395 return 0;
396 else
397 {
398 val -= val & (val-1); // isolate lowest bit
399
400 index_t idx = 0;
401 const index_t nbits = HI - LO;
402
403 if (nbits > 8)
404 {
405 if (val & 0xffffffff00000000ul)
406 idx += 32;
407 if (val & 0xffff0000ffff0000ul)
408 idx += 16;
409 if (val & 0xff00ff00ff00ff00ul)
410 idx += 8;
411 }
412 if (val & 0xf0f0f0f0f0f0f0f0ul)
413 idx += 4;
414 if (val & 0xccccccccccccccccul)
415 idx += 2;
416 if (val & 0xaaaaaaaaaaaaaaaaul)
417 idx += 1;
418
419 return idx + ((idx < -LO) ? LO : LO+1);
420 }
421 }
422#elif (_GLUCAT_BITS_PER_ULONG == 32)
424 template<const index_t LO, const index_t HI>
425 inline
426 index_t
428 min() const
429 {
430 // Reference: [JA], 1.3
431 unsigned long val = bitset_t::to_ulong();
432 if (val == 0)
433 return 0;
434 else
435 {
436 val -= val & (val-1); // isolate lowest bit
437
438 index_t idx = 0;
439 const index_t nbits = HI - LO;
440 if (nbits > 8)
441 {
442 if (val & 0xffff0000ul)
443 idx += 16;
444 if (val & 0xff00ff00ul)
445 idx += 8;
446 }
447 if (val & 0xf0f0f0f0ul)
448 idx += 4;
449 if (val & 0xccccccccul)
450 idx += 2;
451 if (val & 0xaaaaaaaaul)
452 idx += 1;
453
454 return idx + ((idx < -LO) ? LO : LO+1);
455 }
456 }
457#else
459 template<const index_t LO, const index_t HI>
460 auto
462 min() const -> index_t
463 {
464 for (auto
465 idx = LO;
466 idx != 0;
467 ++idx)
468 if (this->test(idx))
469 return idx;
470 for (auto
471 idx = index_t(1);
472 idx <= HI;
473 ++idx)
474 if (this->test(idx))
475 return idx;
476 return 0;
477 }
478#endif
479
480#if (_GLUCAT_BITS_PER_ULONG == 64)
482 template<const index_t LO, const index_t HI>
483 inline
484 auto
486 max() const -> index_t
487 {
488 // Reference: [JA], 1.6
489 auto val = bitset_t::to_ulong();
490 if (val == 0)
491 return 0;
492 else
493 {
494 auto idx = index_t(0);
495 const auto nbits = HI - LO;
496 if (nbits > 8)
497 {
498 if (val & 0xffffffff00000000ul)
499 { val >>= 32; idx += 32; }
500 if (val & 0x00000000ffff0000ul)
501 { val >>= 16; idx += 16; }
502 if (val & 0x000000000000ff00ul)
503 { val >>= 8; idx += 8; }
504 }
505 if (val & 0x00000000000000f0ul)
506 { val >>= 4; idx += 4; }
507 if (val & 0x000000000000000cul)
508 { val >>= 2; idx += 2; }
509 if (val & 0x0000000000000002ul)
510 { idx += 1; }
511 return idx + ((idx < -LO) ? LO : LO+1);
512 }
513 }
514#elif (_GLUCAT_BITS_PER_ULONG == 32)
516 template<const index_t LO, const index_t HI>
517 inline
518 auto
520 max() const -> index_t
521 {
522 // Reference: [JA], 1.6
523 auto val = bitset_t::to_ulong();
524 if (val == 0)
525 return 0;
526 else
527 {
528 auto idx = index_t(0);
529 const auto nbits = HI - LO;
530 if (nbits > 8)
531 {
532 if (val & 0xffff0000ul)
533 { val >>= 16; idx += 16; }
534 if (val & 0x0000ff00ul)
535 { val >>= 8; idx += 8; }
536 }
537 if (val & 0x000000f0ul)
538 { val >>= 4; idx += 4; }
539 if (val & 0x0000000cul)
540 { val >>= 2; idx += 2; }
541 if (val & 0x00000002ul)
542 { idx += 1; }
543 return idx + ((idx < -LO) ? LO : LO+1);
544 }
545 }
546#else
548 template<const index_t LO, const index_t HI>
549 auto
551 max() const -> index_t
552 {
553 for (auto
554 idx = HI;
555 idx != 0;
556 --idx)
557 if (this->test(idx))
558 return idx;
559 for (auto
560 idx = index_t(-1);
561 idx >= LO;
562 --idx)
563 if (this->test(idx))
564 return idx;
565 return 0;
566 }
567#endif
568
570 // eg. {3,4,5} is less than {3,7,8}
571 template<const index_t LO, const index_t HI>
572 inline
573 auto
574 compare(const index_set<LO,HI>& a, const index_set<LO,HI>& b) -> int
575 {
576 return (a == b)
577 ? 0
578 : a.lex_less_than(b)
579 ? -1
580 : 1;
581 }
582
584 // eg. {3,4,5} is less than {3,7,8}
585 template<const index_t LO, const index_t HI>
586 inline
587 auto
589 lex_less_than(const index_set_t rhs) const -> bool
590 { return bitset_t::to_ulong() < rhs.bitset_t::to_ulong(); }
591
593 // Order by count, then order lexicographically within the equivalence class of count.
594 template<const index_t LO, const index_t HI>
595 inline
596 auto
598 operator< (const index_set_t rhs) const -> bool
599 {
600 const auto this_grade = this->count();
601 const auto rhs_grade = rhs.count();
602 return (this_grade < rhs_grade)
603 ? true
604 : (this_grade > rhs_grade)
605 ? false
606 : this->lex_less_than(rhs);
607 }
608
610 template<const index_t LO, const index_t HI>
611 auto
612 operator<< (std::ostream& os, const index_set<LO,HI>& ist) -> std::ostream&
613 {
614 index_t i;
615 os << '{';
616 for (i = LO;
617 (i <= HI) && !(ist[i]);
618 ++i)
619 { }
620 if (i <= HI)
621 os << i;
622 for (++i;
623 i <= HI;
624 ++i)
625 if (ist[i])
626 os << ',' << i;
627 os << '}';
628 return os;
629 }
630
632 template<const index_t LO, const index_t HI>
633 auto
634 operator>> (std::istream& s, index_set<LO,HI>& ist) -> std::istream&
635 {
636 // Parsing variables.
637 auto i = index_t(0);
638 using index_set_t = index_set<LO,HI>;
639 auto local_ist = index_set_t();
640 // Parsing control variables.
641 auto parse_index_list = true;
642 auto expect_closing_brace = false;
643 auto expect_index = false;
644 // Parse an optional opening brace.
645 auto c = s.peek();
646 // If there is a failure or end of file, this ends parsing.
647 if (!s.good())
648 parse_index_list = false;
649 else
650 { // Check for an opening brace.
651 expect_closing_brace = (c == int('{'));
652 if (expect_closing_brace)
653 { // Consume the opening brace.
654 s.get();
655 // The next character may be a closing brace,
656 // indicating the empty index set.
657 c = s.peek();
658 if (s.good() && (c == int('}')))
659 { // A closing brace has been parsed and is no longer expected.
660 expect_closing_brace = false;
661 // Consume the closing brace.
662 s.get();
663 // This ends parsing.
664 parse_index_list = false;
665 }
666 }
667 }
668 if (s.good() && parse_index_list)
669 { // Parse an optional index list.
670 // The index list starts with a first index.
671 for (s >> i;
672 !s.fail();
673 s >> i)
674 { // An index has been parsed. Check to see if it is in range.
675 if ((i < LO) || (i > HI))
676 { // An index out of range is a failure.
677 s.clear(std::istream::failbit);
678 break;
679 }
680 // Add the index to the index set local_ist.
681 local_ist.set(i);
682 // Immediately after parsing an index, an index is no longer expected.
683 expect_index = false;
684 // Reading the index may have resulted in an end of file condition.
685 // If so, this ends the index list.
686 if (s.eof())
687 break;
688 // The index list continues with a comma, and
689 // may be ended by a closing brace, if it was begun with an opening brace.
690 // Parse a possible comma or closing brace.
691 c = s.peek();
692 if (!s.good())
693 break;
694 // First, test for a closing brace, if expected.
695 if (expect_closing_brace && (c == int('}')))
696 { // Consume the closing brace.
697 s.get();
698 // Immediately after parsing the closing brace, it is no longer expected.
699 expect_closing_brace = false;
700 // A closing brace ends the index list.
701 break;
702 }
703 // Now test for a comma.
704 if (c == int(','))
705 { // Consume the comma.
706 s.get();
707 // A index is expected after the comma.
708 expect_index = true;
709 }
710 else
711 { // Any other character here is a failure.
712 s.clear(std::istream::failbit);
713 break;
714 }
715 }
716 }
717 // If an index or a closing brace is still expected, this is a failure.
718 if (expect_index || expect_closing_brace)
719 s.clear(std::istream::failbit);
720 // End of file is not a failure.
721 if (s)
722 { // The index set has been successfully parsed.
723 ist = local_ist;
724 }
725 return s;
726 }
727
729 template<const index_t LO, const index_t HI>
730 inline
731 auto
733 is_contiguous () const -> bool
734 {
735 const auto min_index = this->min();
736 const auto max_index = this->max();
737 return (min_index < 0 && max_index > 0)
738 ? max_index - min_index == this->count()
739 : (min_index == 1 || max_index == -1) &&
740 (max_index - min_index == this->count() - 1);
741 }
742
744 template<const index_t LO, const index_t HI>
745 inline
746 auto
748 fold() const -> const
749 index_set<LO,HI>
750 { return this->fold(*this, true); }
751
753 template<const index_t LO, const index_t HI>
754 auto
756 fold(const index_set_t frm, const bool prechecked) const -> const
758 {
759 if (!prechecked && ((*this | frm) != frm))
760 throw error_t("fold(frm): cannot fold from outside of frame");
761 const auto frm_min = frm.min();
762 const auto frm_max = frm.max();
763 auto result = index_set_t();
764 auto fold_idx = index_t(-1);
765 for (auto
766 unfold_idx = fold_idx;
767 unfold_idx >= frm_min;
768 --unfold_idx)
769 if (frm.test(unfold_idx))
770 // result.set(fold_idx--, this->test(unfold_idx));
771 {
772 if (this->test(unfold_idx))
773 result.set(fold_idx);
774 --fold_idx;
775 }
776 fold_idx = index_t(1);
777 for (auto
778 unfold_idx = fold_idx;
779 unfold_idx <= frm_max;
780 ++unfold_idx)
781 if (frm.test(unfold_idx))
782 // result.set(fold_idx++, this->test(unfold_idx));
783 {
784 if (this->test(unfold_idx))
785 result.set(fold_idx);
786 ++fold_idx;
787 }
788 return result;
789 }
790
792 template<const index_t LO, const index_t HI>
793 auto
795 unfold(const index_set_t frm, const bool prechecked) const -> const index_set_t
796 {
797 const char* msg =
798 "unfold(frm): cannot unfold into a smaller frame";
799 const auto frm_min = frm.min();
800 const auto frm_max = frm.max();
801 auto result = index_set_t();
802 auto fold_idx = index_t(-1);
803 for (auto
804 unfold_idx = fold_idx;
805 unfold_idx >= frm_min;
806 --unfold_idx)
807 if (frm.test(unfold_idx))
808 if (this->test(fold_idx--))
809 result.set(unfold_idx);
810 if (!prechecked && ((fold_idx+1) > this->min()))
811 throw error_t(msg);
812 fold_idx = index_t(1);
813 for (auto
814 unfold_idx = fold_idx;
815 unfold_idx <= frm_max;
816 ++unfold_idx)
817 if (frm.test(unfold_idx))
818 if (this->test(fold_idx++))
819 result.set(unfold_idx);
820 if (!prechecked && ((fold_idx-1) < this->max()))
821 throw error_t(msg);
822 return result;
823 }
824
826 template<const index_t LO, const index_t HI>
827 inline
828 auto
830 value_of_fold(const index_set_t frm) const -> set_value_t
831 {
832 const auto min_index = frm.fold().min();
833 if (min_index == 0)
834 return 0;
835 else
836 {
837 const auto folded_set = this->fold(frm);
838 const auto skip = min_index > 0 ? index_t(1) : index_t(0);
839 return folded_set.bitset_t::to_ulong() >> (min_index-LO-skip);
840 }
841 }
842
844 inline
845 static
846 auto inverse_reversed_gray(unsigned long x) -> unsigned long
847 {
848 // Reference: [JA]
849#if (_GLUCAT_BITS_PER_ULONG >= 64)
850 x ^= x << 32; // for 64-bit words
851#endif
852 x ^= x << 16; // reversed_gray ** 16
853 x ^= x << 8; // reversed_gray ** 8
854 x ^= x << 4; // reversed_gray ** 4
855 x ^= x << 2; // reversed_gray ** 2
856 x ^= x << 1; // reversed_gray ** 1
857 return x;
858 }
859
861 inline
862 static
863 auto inverse_gray(unsigned long x) -> unsigned long
864 {
865 // Reference: [JA]
866#if (_GLUCAT_BITS_PER_ULONG >= 64)
867 x ^= x >> 32; // for 64-bit words
868#endif
869 x ^= x >> 16; // gray ** 16
870 x ^= x >> 8; // gray ** 8
871 x ^= x >> 4; // gray ** 4
872 x ^= x >> 2; // gray ** 2
873 x ^= x >> 1; // gray ** 1
874 return x;
875 }
876
878 template<const index_t LO, const index_t HI>
879 auto
881 sign_of_mult(const index_set_t rhs) const -> int
882 {
883 // Implemented using Walsh functions and Gray codes.
884 // Reference: [L] Chapter 21, 21.3
885 // Reference: [JA]
886 const auto uthis = this->bitset_t::to_ulong();
887 const auto urhs = rhs.bitset_t::to_ulong();
888 const auto nbits = HI - LO;
889 auto negative = 0UL;
890 if (nbits > 8)
891 {
892 // Set h to be the inverse reversed Gray code of rhs.
893 // This sets each bit of h to be the cumulative ^ of
894 // the same and lower bits of rhs.
895 const auto h = inverse_reversed_gray(urhs);
896 // Set k to be the inverse Gray code of *this & h.
897 // This sets the low bit of k to be parity(*this & h).
898 const auto k = inverse_gray(uthis & h);
899 // Set q to be the inverse Gray code of the positive part of *this & rhs.
900 const auto q = inverse_gray((uthis & urhs) >> -LO);
901 negative = k ^ q;
902 }
903 else
904 {
905 auto h = 0UL;
906 for (auto
907 j = index_t(0);
908 j < -LO;
909 ++j)
910 {
911 h ^= urhs >> j;
912 negative ^= h & (uthis >> j);
913 }
914 for (auto
915 j = index_t(-LO);
916 j < nbits;
917 ++j)
918 {
919 negative ^= h & (uthis >> j);
920 h ^= urhs >> j;
921 }
922 }
923 return 1 - int((negative & 1) << 1);
924 }
925
927 template<const index_t LO, const index_t HI>
928 inline
929 auto
931 sign_of_square() const -> int
932 {
933 auto result = 1 - int((this->count_neg() % 2) << 1);
934 switch (this->count() % 4)
935 {
936 case 2:
937 case 3:
938 result *= -1;
939 break;
940 default:
941 break;
942 }
943 return result;
944 }
945
947 template<const index_t LO, const index_t HI>
948 inline
949 auto
951 hash_fn() const -> size_t
952 {
953 static const auto lo_mask = (1UL << -LO) - 1UL;
954 const auto uthis = bitset_t::to_ulong();
955 const auto neg_part = uthis & lo_mask;
956 const auto pos_part = uthis >> -LO;
957 return size_t(neg_part ^ pos_part);
958 }
959
961 inline
962 auto
964 { return (j < 0) ? -1 : 1; }
965
967 template<const index_t LO, const index_t HI>
968 inline
969 auto
971 { return std::min(ist.min(), 0); }
972
974 template<const index_t LO, const index_t HI>
975 inline
976 auto
978 { return std::max(ist.max(), 0); }
979
980// index_set reference
981
983 template<const index_t LO, const index_t HI>
984 inline
986 reference( index_set_t& ist, index_t idx ) :
987 m_pst(&ist),
988 m_idx(idx)
989 { }
990
992 template<const index_t LO, const index_t HI>
993 inline
994 auto
996 operator== (const reference& c_j) const -> bool
997 { return m_pst == c_j.m_pst && m_idx == c_j.m_idx; }
998
1000 template<const index_t LO, const index_t HI>
1001 inline
1002 auto
1004 operator= (bool x) -> reference&
1005 {
1006 if ( x )
1007 m_pst->set(m_idx);
1008 else
1009 m_pst->reset(m_idx);
1010 return *this;
1011 }
1012
1014 template<const index_t LO, const index_t HI>
1015 inline
1016 auto
1018 operator= (const reference& c_j) -> reference&
1019 {
1020 if (&c_j != this && c_j != *this)
1021 {
1022 if ( (c_j.m_pst)[c_j.m_idx] )
1023 m_pst->set(m_idx);
1024 else
1025 m_pst->reset(m_idx);
1026 }
1027 return *this;
1028 }
1029
1031 template<const index_t LO, const index_t HI>
1032 inline
1033 auto
1035 operator~ () const -> bool
1036 { return !(m_pst->test(m_idx)); }
1037
1039 template<const index_t LO, const index_t HI>
1040 inline
1042 operator bool () const
1043 { return m_pst->test(m_idx); }
1044
1046 template<const index_t LO, const index_t HI>
1047 inline
1048 auto
1050 flip() -> reference&
1051 {
1052 m_pst->flip(m_idx);
1053 return *this;
1054 }
1055}
1056#endif // _GLUCAT_INDEX_SET_IMP_H
Index set member reference.
Definition index_set.h:177
auto operator=(const bool x) -> reference &
for b[i] = x;
auto flip() -> reference &
for b[i].flip();
reference()=delete
Default constructor is deleted.
auto operator~() const -> bool
Flips a bit.
auto operator==(const reference &c_j) const -> bool
for b[i] == c[j];
Index set class based on std::bitset<> in Gnu standard C++ library.
Definition index_set.h:75
auto value_of_fold(const index_set_t frm) const -> set_value_t
The set value of the fold of this index set within the given frame.
auto count() const -> index_t
Cardinality: Number of indices included in set.
std::bitset< HI - LO > bitset_t
Definition index_set.h:81
auto flip() -> index_set_t &
Set complement, except 0: flip all bits, except 0.
auto operator[](const index_t idx) const -> bool
Subscripting: Test idx for membership: test value of bit idx.
error< index_set > error_t
Definition index_set.h:82
auto lex_less_than(const index_set_t rhs) const -> bool
Lexicographic ordering of two sets: *this < rhs.
auto reset() -> index_set_t &
Make set empty: Set all bits to 0.
auto sign_of_square() const -> int
Sign of geometric square of a Clifford basis element.
auto operator<(const index_set_t rhs) const -> bool
Less than operator used for comparisons, map, etc.
auto min() const -> index_t
Minimum member.
auto hash_fn() const -> size_t
Hash function.
auto set() -> index_set_t &
Include all indices except 0: set all bits except 0.
auto is_contiguous() const -> bool
Determine if the index set is contiguous, ie. has no gaps.
auto count_pos() const -> index_t
Number of positive indices included in set.
std::pair< index_t, index_t > index_pair_t
Definition index_set.h:85
auto operator|=(const index_set_t rhs) -> index_set_t &
Set union: or.
auto count_neg() const -> index_t
Number of negative indices included in set.
auto operator!=(const index_set_t rhs) const -> bool
Inequality.
static auto classname() -> const std::string
auto operator~() const -> index_set_t
Set complement: not.
index_set index_set_t
Definition index_set.h:84
auto operator&=(const index_set_t rhs) -> index_set_t &
Set intersection: and.
auto operator^=(const index_set_t rhs) -> index_set_t &
Symmetric set difference: exclusive or.
auto max() const -> index_t
Maximum member.
auto sign_of_mult(const index_set_t ist) const -> int
Sign of geometric product of two Clifford basis elements.
auto unfold(const index_set_t frm, const bool prechecked=false) const -> const index_set_t
Unfold this index set within the given frame.
index_set()=default
Default constructor creates an empty set.
auto operator==(const index_set_t rhs) const -> bool
Equality.
auto test(const index_t idx) const -> bool
Test idx for membership: test value of bit idx.
auto fold() const -> const index_set_t
Fold this index set within itself as a frame.
auto operator<<(std::ostream &os, const framed_multi< Scalar_T, LO, HI, Tune_P > &val) -> std::ostream &
Write multivector to output.
auto operator|(const Multivector< Scalar_T, LO, HI, Tune_P > &lhs, const RHS< Scalar_T, LO, HI, Tune_P > &rhs) -> const Multivector< Scalar_T, LO, HI, Tune_P >
Transformation via twisted adjoint action.
auto operator&(const Multivector< Scalar_T, LO, HI, Tune_P > &lhs, const RHS< Scalar_T, LO, HI, Tune_P > &rhs) -> const Multivector< Scalar_T, LO, HI, Tune_P >
Inner product.
auto compare(const index_set< LO, HI > &a, const index_set< LO, HI > &b) -> int
"lexicographic compare" eg. {3,4,5} is less than {3,7,8}
static auto inverse_reversed_gray(unsigned long x) -> unsigned long
Inverse reversed Gray code.
static auto inverse_gray(unsigned long x) -> unsigned long
Inverse Gray code.
auto sign_of_square(index_t j) -> int
Square of generator {j}.
auto min_neg(const index_set< LO, HI > &ist) -> index_t
Minimum negative index, or 0 if none.
unsigned long set_value_t
Size of set_value_t should be enough to contain index_set<LO,HI>
Definition global.h:79
auto operator>>(std::istream &s, framed_multi< Scalar_T, LO, HI, Tune_P > &val) -> std::istream &
Read multivector from input.
int index_t
Size of index_t should be enough to represent LO, HI.
Definition global.h:77
auto operator^(const Multivector< Scalar_T, LO, HI, Tune_P > &lhs, const RHS< Scalar_T, LO, HI, Tune_P > &rhs) -> const Multivector< Scalar_T, LO, HI, Tune_P >
Outer product.
auto max_pos(const index_set< LO, HI > &ist) -> index_t
Maximum positive index, or 0 if none.