Frobby 0.9.5
IdealTree.cpp
Go to the documentation of this file.
1/* Frobby: Software for monomial ideal computations.
2 Copyright (C) 2007 Bjarke Hammersholt Roune (www.broune.com)
3 Copyright (C) 2010 University of Aarhus
4 Contact Bjarke Hammersholt Roune for license information (www.broune.com)
5
6 This program is free software; you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation; either version 2 of the License, or
9 (at your option) any later version.
10
11 This program is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with this program. If not, see http://www.gnu.org/licenses/.
18*/
19#include "stdinc.h"
20#include "IdealTree.h"
21
22#include "Ideal.h"
23#include "Term.h"
24
25namespace {
26 const size_t MaxLeafSize = 60;
27
30 const Exponent* term,
31 size_t varCount) {
32 for (; begin != end; ++begin)
33 if (Term::strictlyDivides(*begin, term, varCount))
34 return true;
35 return false;
36 }
37}
38
40public:
43 size_t varCount):
44 _begin(begin), _end(end), _varCount(varCount) {
45 makeTree();
46 }
47
48 void makeTree();
49 bool strictlyContains(const Exponent* term) const;
50 size_t getVarCount() const {return _varCount;}
51
52private:
57 size_t _varCount;
58 size_t _var;
59 size_t _pivot;
60};
61
63 if ((size_t)distance(_begin, _end) <= MaxLeafSize)
64 return;
65
66 Term lcm(_varCount);
67 Term gcd(*_begin, _varCount);
68 for (Ideal::const_iterator it = _begin; it != _end; ++it) {
69 lcm.lcm(lcm, *it);
70 gcd.gcd(gcd, *it);
71 }
72
73 while (true) {
74 size_t maxVar = 0;
75 for (size_t var = 1; var < _varCount; ++var)
76 if (lcm[var] - gcd[var] > lcm[maxVar] - gcd[maxVar])
77 maxVar = var;
78
79 // TODO: could this ever happen?
80 if (lcm[maxVar] == 0)
81 break; // we are not making any progress anyway.
82
83 ASSERT(lcm[maxVar] >= 1);
84 _var = maxVar;
85
86 // It's significant that we are rounding down to ensure that
87 // neither of _lessOrEqual and _greater becomes empty.
88 _pivot = (lcm[maxVar] + gcd[maxVar]) >> 1; // Note: x >> 1 == x / 2
89
92
93 // put those strictly divisible by _var^_pivot to right and the
94 // rest to the left.
95 while (left != right) {
96 // Find left-most element that should go right.
97 while ((*left)[_var] <= _pivot && left != right)
98 ++left;
99
100 // Find right-most element that should go left.
101 while ((*right)[_var] > _pivot && left != right)
102 --right;
103
104 // Swap the two found elements so that they both go into the
105 // right place.
106 using std::swap;
107 swap(*left, *right);
108 }
109 ASSERT((*(_end - 1))[_var] > _pivot);
110 ASSERT((*_begin)[_var] <= _pivot);
111
112 // Make middle the first element that went right, so that the two
113 // ranges become [_begin, middle) and [middle, _end).
114 ASSERT((*right)[_var] > _pivot)
116 while ((*middle)[_var] > _pivot) {
117 ASSERT(middle != _begin);
118 --middle;
119 }
120 ++middle;
121 ASSERT(_begin < middle && middle <= _end);
122 ASSERT((*middle)[_var] > _pivot);
123 ASSERT((*(middle - 1))[_var] <= _pivot);
124
126 _greater.reset(new Node(middle, _end, _varCount));
127
128 _lessOrEqual->makeTree();
129 _greater->makeTree();
130 return;
131 }
132}
133
135 if (_lessOrEqual.get() != 0) {
136 ASSERT(_greater.get() != 0);
137 bool returnValue =
138 _lessOrEqual->strictlyContains(term) ||
139 _greater->strictlyContains(term);
140 ASSERT(returnValue == rawStrictlyDivides(_begin, _end, term, _varCount));
141 return returnValue;
142 } else {
143 ASSERT(_greater.get() == 0);
144 return rawStrictlyDivides(_begin, _end, term, _varCount);
145 }
146}
147
149 // not using initialization for this to avoid depending on order of
150 // initialization of members.
151 _storage.reset(new Ideal(ideal));
152 _root.reset
153 (new Node(_storage->begin(), _storage->end(), ideal.getVarCount()));
154}
155
157 // Destructor defined so auto_ptr<T> in the header does not need
158 // definition of T.
159}
160
161bool IdealTree::strictlyContains(const Exponent* term) const {
162 ASSERT(_root.get() != 0);
163 return _root->strictlyContains(term);
164}
165
167 ASSERT(_root.get() != 0);
168 return _root->getVarCount();
169}
void swap(HilbertSlice &a, HilbertSlice &b)
void nameFactoryRegister(NameFactory< AbstractProduct > &factory)
Registers the string returned by ConcreteProduct::getStaticName() to a function that default-construc...
bool strictlyContains(const Exponent *term) const
Ideal::iterator _begin
Definition IdealTree.cpp:55
size_t getVarCount() const
Definition IdealTree.cpp:50
auto_ptr< Node > _lessOrEqual
Definition IdealTree.cpp:53
Ideal::iterator _end
Definition IdealTree.cpp:56
auto_ptr< Node > _greater
Definition IdealTree.cpp:54
Node(Ideal::iterator begin, Ideal::iterator end, size_t varCount)
Definition IdealTree.cpp:41
IdealTree(const Ideal &ideal)
auto_ptr< Ideal > _storage
Definition IdealTree.h:40
bool strictlyContains(const Exponent *term) const
size_t getVarCount() const
auto_ptr< Node > _root
Definition IdealTree.h:41
Represents a monomial ideal with int exponents.
Definition Ideal.h:27
Cont::const_iterator const_iterator
Definition Ideal.h:43
Cont::iterator iterator
Definition Ideal.h:44
size_t getVarCount() const
Definition Ideal.h:56
Term represents a product of variables which does not include a coefficient.
Definition Term.h:49
static bool strictlyDivides(const Exponent *a, const Exponent *b, size_t varCount)
Returns whether a strictly divides b.
Definition Term.h:196
unsigned int Exponent
Definition stdinc.h:89
#define ASSERT(X)
Definition stdinc.h:86