Frobby 0.9.5
SliceFacade.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
4 This program is free software; you can redistribute it and/or modify
5 it under the terms of the GNU General Public License as published by
6 the Free Software Foundation; either version 2 of the License, or
7 (at your option) any later version.
8
9 This program is distributed in the hope that it will be useful,
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 GNU General Public License for more details.
13
14 You should have received a copy of the GNU General Public License
15 along with this program. If not, see http://www.gnu.org/licenses/.
16*/
17#include "stdinc.h"
18#include "SliceFacade.h"
19
20#include "BigTermConsumer.h"
21#include "CoefBigTermConsumer.h"
22#include "TermTranslator.h"
23#include "BigIdeal.h"
24#include "Ideal.h"
25#include "Term.h"
26#include "MsmStrategy.h"
29#include "DebugStrategy.h"
30#include "DecomRecorder.h"
31#include "TermGrader.h"
32#include "OptimizeStrategy.h"
34#include "HilbertStrategy.h"
35#include "IOHandler.h"
36#include "BigPolynomial.h"
38#include "CoefBigTermRecorder.h"
40#include "VarSorter.h"
41#include "StatisticsStrategy.h"
43#include "SizeMaxIndepSetAlg.h"
44#include "SliceParams.h"
45#include "error.h"
46#include "display.h"
47
48#include <iterator>
49
50SliceFacade::SliceFacade(const SliceParams& params, const DataType& output):
51 Facade(params.getPrintActions()),
52 _params(params) {
54 _common.readIdealAndSetOutput(params, output);
55}
56
58 const BigIdeal& ideal,
59 BigTermConsumer& consumer):
60 Facade(params.getPrintActions()),
61 _params(params) {
63 _common.setIdealAndIdealOutput(params, ideal, consumer);
64}
65
67 const BigIdeal& ideal,
68 CoefBigTermConsumer& consumer):
69 Facade(params.getPrintActions()),
70 _params(params) {
72 _common.setIdealAndPolyOutput(params, ideal, consumer);
73}
74
76}
77
80 beginAction("Computing multigraded Hilbert-Poincare series.");
81
82 auto_ptr<CoefTermConsumer> consumer = _common.makeTranslatedPolyConsumer();
83
84 consumer->consumeRing(_common.getNames());
85 consumer->beginConsuming();
86 HilbertStrategy strategy(consumer.get(), _split.get());
88 consumer->doneConsuming();
89
90 endAction();
91}
92
95 beginAction("Computing univariate Hilbert-Poincare series.");
96
97 auto_ptr<CoefTermConsumer> consumer =
99
100 consumer->consumeRing(_common.getNames());
101 consumer->beginConsuming();
102 HilbertStrategy strategy(consumer.get(), _split.get());
104 consumer->doneConsuming();
105
106 endAction();
107}
108
112}
113
116
118 if (codimension) {
119 // convert to mpz_class before increment to ensure no overflow.
120 return mpz_class(_common.getIdeal().getVarCount()) + 1;
121 } else
122 return -1;
123 }
124
125 // todo: inline?
126 takeRadical();
127
128 beginAction("Preparing to compute dimension.");
129
130 vector<mpz_class> v;
131 fill_n(back_inserter(v), _common.getIdeal().getVarCount(), -1);
132
133 endAction();
134
135 mpz_class minusCodimension;
136#ifdef DEBUG
137 // Only define hasComponents when DEBUG is defined since otherwise
138 // GCC will warn about hasComponents not being used.
139 bool hasComponents =
140#endif
141 solveIrreducibleDecompositionProgram(v, minusCodimension, false);
142 ASSERT(hasComponents);
143
144 if (codimension)
145 return -minusCodimension;
146 else
147 return v.size() + minusCodimension;
148}
149
152
153 size_t varCount = _common.getIdeal().getVarCount();
154
155 Ideal irreducibleDecom(varCount);
156 {
157 DecomRecorder recorder(&irreducibleDecom);
158 produceEncodedIrrDecom(recorder);
159 }
160
162 ("Computing primary decomposition from irreducible decomposition.");
163
164 // Do intersection of each component also using irreducible
165 // decomposition of the dual. We can't use the Alexander dual
166 // methods, since those switch around the translator to emit altered
167 // big integers, while keeping the small integers the same, but we
168 // want to keep this in small integers. So we have to do the dual
169 // thing here.
170
171 // To get actual supports.
173
174 // To collect same-support vectors together.
175 irreducibleDecom.sortReverseLex();
176
177 Term lcm(varCount);
178 irreducibleDecom.getLcm(lcm);
179
180 Term tmp(varCount);
181 Term support(varCount);
182
184 Ideal& primaryComponentDual = _common.getIdeal();
185 Ideal primaryComponent(varCount);
186
187 DecomRecorder recorder(&primaryComponent);
188
189 auto_ptr<TermConsumer> consumer = _common.makeTranslatedIdealConsumer();
190 consumer->consumeRing(_common.getNames());
191 consumer->beginConsumingList();
192
193 Ideal::const_iterator stop = irreducibleDecom.end();
194 Ideal::const_iterator it = irreducibleDecom.begin();
195 while (it != stop) {
196 // Get all vectors with same support.
197 support = *it;
198 do {
199 tmp.encodedDual(*it, lcm);
200 primaryComponentDual.insert(tmp);
201 ++it;
202 } while (it != stop && support.hasSameSupport(*it));
203 ASSERT(!primaryComponentDual.isZeroIdeal());
204
205 _common.getTranslator().addPurePowersAtInfinity(primaryComponentDual);
206 {
207 MsmStrategy strategy(&recorder, _split.get());
209 }
211
212 consumer->beginConsuming();
213 for (Ideal::const_iterator dualTerm = primaryComponent.begin();
214 dualTerm != primaryComponent.end(); ++dualTerm) {
215 tmp.encodedDual(*dualTerm, lcm);
216 consumer->consume(tmp);
217 }
218 consumer->doneConsuming();
219
220 primaryComponent.clear();
221 primaryComponentDual.clear();
222 }
223
224 consumer->doneConsumingList();
225
226 endAction();
227}
228
231 beginAction("Computing maximal staircase monomials.");
232
233 auto_ptr<TermConsumer> consumer = _common.makeTranslatedIdealConsumer();
234 consumer->consumeRing(_common.getNames());
235 MsmStrategy strategy(consumer.get(), _split.get());
237
238 endAction();
239}
240
243
244 beginAction("Preparing to compute maximal standard monomials.");
246 endAction();
248}
249
250void SliceFacade::computeAlexanderDual(const vector<mpz_class>& point) {
252 ASSERT(point.size() == _common.getIdeal().getVarCount());
253
254 beginAction("Ensuring specified point is divisible by lcm.");
255 vector<mpz_class> lcm;
257
258 for (size_t var = 0; var < lcm.size(); ++var) {
259 if (lcm[var] > point[var]) {
260 endAction();
262 ("The specified point to dualize on is not divisible by the "
263 "least common multiple of the minimal generators of the ideal.");
264 }
265 }
266 endAction();
267
268 beginAction("Preparing to compute Alexander dual.");
270 endAction();
271
273}
274
277
278 beginAction("Computing lcm for Alexander dual.");
279 vector<mpz_class> lcm;
281 endAction();
282
284}
285
288
289 size_t varCount = _common.getIdeal().getVarCount();
290
291 // Obtain generators of radical from irreducible decomposition.
292 Ideal radical(varCount);
293 {
294 Ideal decom(varCount);
295 DecomRecorder recorder(&decom);
296 produceEncodedIrrDecom(recorder);
297
298 beginAction("Computing associated primes from irreducible decomposition.");
299
300 Term tmp(varCount);
301 Ideal::const_iterator stop = decom.end();
302 for (Ideal::const_iterator it = decom.begin(); it != stop; ++it) {
303 for (size_t var = 0; var < varCount; ++var) {
304 // We cannot just check whether (*it)[var] == 0, since the
305 // added fake pure powers map to zero but are not themselves
306 // zero.
307 if (_common.getTranslator().getExponent(var, (*it)[var]) == 0)
308 tmp[var] = 0;
309 else
310 tmp[var] = 1;
311 }
312 radical.insert(tmp);
313 }
314 }
315
316 radical.removeDuplicates();
317
318
319 // Output associated primes.
321 auto_ptr<TermConsumer> consumer = _common.makeTranslatedIdealConsumer();
322
323 consumer->consumeRing(_common.getNames());
324 consumer->beginConsuming();
325 Term tmp(varCount);
326 Ideal::const_iterator stop = radical.end();
327 for (Ideal::const_iterator it = radical.begin(); it != stop; ++it) {
328 tmp = *it;
329 consumer->consume(tmp);
330 }
331 consumer->doneConsuming();
332
333 endAction();
334}
335
337(const vector<mpz_class>& grading,
338 mpz_class& optimalValue,
339 bool reportAllSolutions) {
341 ASSERT(grading.size() == _common.getIdeal().getVarCount());
342
343 beginAction("Preparing to solve optimization program.");
346 endAction();
347
348 return solveProgram(grading, optimalValue, reportAllSolutions);
349}
350
352(const vector<mpz_class>& grading,
353 mpz_class& optimalValue,
354 bool reportAllSolutions) {
356 ASSERT(grading.size() == _common.getIdeal().getVarCount());
357
359 return solveProgram(grading, optimalValue, reportAllSolutions);
360}
361
364 beginAction("Computing irreducible decomposition.");
365
367 MsmStrategy strategy(&consumer, _split.get());
368
369 consumer.consumeRing(_common.getNames());
371
372 endAction();
373}
374
375bool SliceFacade::solveProgram(const vector<mpz_class>& grading,
376 mpz_class& optimalValue,
377 bool reportAllSolutions) {
379 ASSERT(grading.size() == _common.getIdeal().getVarCount());
380
383 ("Turning off Independence splits as they are not supported\n"
384 "for optimization.");
386 }
387
391 ("Bound simplification requires using the bound to eliminate\n"
392 "non-improving slices, which has been turned off. Am now turning\n"
393 "this on.");
395 }
396
397 beginAction("Solving optimization program.");
398
403 } else if (_params.getUseBoundElimination())
405 else
406 boundSetting = OptimizeStrategy::DoNotUseBound;
407
408 TermGrader grader(grading, _common.getTranslator());
409 OptimizeStrategy strategy
410 (grader, _split.get(), reportAllSolutions, boundSetting);
412
413 endAction();
414
415 const Ideal& solution = strategy.getMaximalSolutions();
416
417 auto_ptr<TermConsumer> consumer = _common.makeTranslatedIdealConsumer();
418 consumer->consumeRing(_common.getNames());
419 consumer->consume(solution);
420
421 if (solution.isZeroIdeal())
422 return false;
423 else {
424 optimalValue = strategy.getMaximalValue();
425 return true;
426 }
427}
428
430 return _common.hasIdeal();
431}
432
435
436 beginAction("Taking radical of ideal.");
437
438 bool skip = false;
441 if (lcm.isSquareFree())
442 skip = true;
443
444 if (!skip) {
448 }
449
451
452 endAction();
453}
454
455void SliceFacade::getLcmOfIdeal(vector<mpz_class>& bigLcm) {
457
460
461 bigLcm.clear();
462 bigLcm.reserve(_common.getIdeal().getVarCount());
463 for (size_t var = 0; var < _common.getIdeal().getVarCount(); ++var)
464 bigLcm.push_back(_common.getTranslator().getExponent(var, lcm));
465}
466
471
472 SliceStrategy* strategyWithOptions = &strategy;
473
474 auto_ptr<SliceStrategy> debugStrategy;
475 if (_params.getPrintDebug()) {
476 debugStrategy.reset
477 (new DebugStrategy(strategyWithOptions, stderr));
478 strategyWithOptions = debugStrategy.get();
479 }
480
481 auto_ptr<SliceStrategy> statisticsStrategy;
483 statisticsStrategy.reset
484 (new StatisticsStrategy(strategyWithOptions, stderr));
485 strategyWithOptions = statisticsStrategy.get();
486 }
487
488 ASSERT(strategyWithOptions != 0);
489 strategyWithOptions->run(_common.getIdeal());
490}
void setToZeroOne(TermTranslator &translator)
TermTranslator & getTranslator()
auto_ptr< TermConsumer > makeTranslatedIdealConsumer(bool split=false)
void setIdealAndIdealOutput(const CommonParams &params, const BigIdeal &input, BigTermConsumer &output)
Use given ideal and support ideal output.
auto_ptr< CoefTermConsumer > makeToUnivariatePolyConsumer()
const VarNames & getNames()
void readIdealAndSetOutput(const CommonParams &params, const DataType &output)
Read input ideal and support specified kind of output.
void setIdealAndPolyOutput(const CommonParams &params, const BigIdeal &input, CoefBigTermConsumer &output)
Use given ideal and support polynomial output.
auto_ptr< CoefTermConsumer > makeTranslatedPolyConsumer()
bool getPrintStatistics() const
Returns whether to print statistics on what the algorithm did to standard error after it has run.
Definition: CommonParams.h:66
bool getPrintDebug() const
Returns whether to print information about what the algorithm is doing to standard error as it runs.
Definition: CommonParams.h:61
The intention of this class is to describe the different kinds of mathematical structures that Frobby...
Definition: DataType.h:29
This is the super class of all facades.
Definition: Facade.h:32
void beginAction(const char *message)
Prints message to standard error if printing is turned on, and records the time when the action start...
Definition: Facade.cpp:38
void endAction()
Prints to standard error the time since the last call to beginAction.
Definition: Facade.cpp:51
Represents a monomial ideal with int exponents.
Definition: Ideal.h:27
void removeDuplicates()
Definition: Ideal.cpp:634
void minimize()
Definition: Ideal.cpp:501
bool isZeroIdeal() const
Definition: Ideal.cpp:86
bool containsIdentity() const
Definition: Ideal.cpp:65
void getLcm(Exponent *lcm) const
Sets lcm to the least common multiple of all generators.
Definition: Ideal.cpp:157
void takeRadicalNoMinimize()
Replaces all generators with their support and does not remove any non-minimal generators this may pr...
Definition: Ideal.cpp:660
Cont::const_iterator const_iterator
Definition: Ideal.h:43
void sortReverseLex()
Definition: Ideal.cpp:510
void clear()
Definition: Ideal.cpp:641
void insert(const Exponent *term)
Definition: Ideal.cpp:455
const_iterator end() const
Definition: Ideal.h:49
const_iterator begin() const
Definition: Ideal.h:48
size_t getVarCount() const
Definition: Ideal.h:56
OptimizeStrategy optimizes a function on the maximal standard monomials of a monomial ideal using bra...
BoundSetting
The values of BoundSetting indicate how to use the bound.
@ DoNotUseBound
Make no use of the bound.
@ UseBoundToEliminateAndSimplify
Eliminate non-improving slices and simplify slices by trying to generate non-improving slices that ar...
@ UseBoundToEliminate
Eliminate non-improving slices, achieving a branch-and-bound algorithm in place of the usual backtrac...
const Ideal & getMaximalSolutions()
Returns one of or all of the msm's with optimal value found so far, depending on the value of reportA...
const mpz_class & getMaximalValue()
The optimal value associated to all entries from getMaximalSolutions().
void takeRadical()
void computeAssociatedPrimes()
Compute the associated primes of the ideal.
void runSliceAlgorithmWithOptions(SliceStrategy &strategy)
void computeMultigradedHilbertSeries()
Compute the numerator of the multigraded Hilbert-Poincare series.
Definition: SliceFacade.cpp:78
bool solveStandardProgram(const vector< mpz_class > &grading, mpz_class &value, bool reportAllSolutions)
Solve an optimization program over maximal standard monomials.
bool solveIrreducibleDecompositionProgram(const vector< mpz_class > &grading, mpz_class &optimalValue, bool reportAllSolutions)
Solve an optimization program over irreducible components.
bool isFirstComputation() const
CommonParamsHelper _common
Definition: SliceFacade.h:224
void getLcmOfIdeal(vector< mpz_class > &lcm)
mpz_class computeDimension(bool codimension=false)
Compute the Krull dimension of ideal.
SliceFacade(const SliceParams &params, const DataType &output)
Definition: SliceFacade.cpp:50
void computeUnivariateHilbertSeries()
Compute the numerator of the univariate Hilbert-Poincare series.
Definition: SliceFacade.cpp:93
SliceParams _params
Definition: SliceFacade.h:223
auto_ptr< SplitStrategy > _split
Definition: SliceFacade.h:225
void computePrimaryDecomposition()
Compute the unique "nicest" primary decomposition of the ideal.
bool solveProgram(const vector< mpz_class > &grading, mpz_class &optimalValue, bool reportAllSolutions)
void computeIrreducibleDecomposition(bool encode)
Compute the unique irredundant set of irreducible ideals whose intersection equals ideal.
void computeMaximalStaircaseMonomials()
Compute the maximal staircase monomials of the ideal.
void computeMaximalStandardMonomials()
Compute the maximal standard monomials of the ideal.
void produceEncodedIrrDecom(TermConsumer &consumer)
void computeAlexanderDual()
Compute the Alexander dual of the ideal.
bool getUseSimplification() const
Apply simplification to the state of the algorithm when possible.
void useIndependenceSplits(bool value)
Definition: SliceParams.h:34
const string & getSplit() const
Definition: SliceParams.h:30
bool getUseIndependenceSplits() const
Definition: SliceParams.h:33
void useBoundElimination(bool value)
Definition: SliceParams.h:39
bool getUseBoundSimplification() const
Returns whether to simplify slices by seeking to generate non-improving slices that are then eliminat...
Definition: SliceParams.h:44
bool getUseBoundElimination() const
Returns whether to use branch-and-bound to speed up Slice optimization computations by eliminating no...
Definition: SliceParams.h:38
This class describes the interface of a strategy object for the Slice Algorithm.
Definition: SliceStrategy.h:33
virtual void setUseSimplification(bool use)=0
This method should only be called before calling run().
virtual void setUseIndependence(bool use)=0
This method should only be called before calling run().
virtual void run(const Ideal &ideal)=0
Run the Slice algorithm.
static auto_ptr< SplitStrategy > createStrategy(const string &prefix)
Returns the strategy whose name has the given prefix.
A wrapper for a SliceStrategy that collects statistics on what is going on, while delegating everythi...
This class is used to transfer terms one at a time from one part of the program to another,...
Definition: TermConsumer.h:36
virtual void consumeRing(const VarNames &names)
Tell the consumer which ring is being used.
A TermGrader assigns a value, the degree, to each monomial.
Definition: TermGrader.h:27
const mpz_class & getExponent(size_t variable, Exponent exponent) const
This method translates from IDs to arbitrary precision integers.
void addPurePowersAtInfinity(Ideal &ideal) const
Adds a generator of the form v^e, e > 0, for any variable v where generator of that form is not alrea...
void dualize(const vector< mpz_class > &a)
Replaces var^v by var^(a[i] - v) except that var^0 is left alone.
void setInfinityPowersToZero(Ideal &ideal) const
The method addPurePowersAtInfinity adds high exponents that map to zero.
void decrement()
Replaces var^v by var^(v-1).
Term represents a product of variables which does not include a coefficient.
Definition: Term.h:49
static bool hasSameSupport(const Exponent *a, const Exponent *b, size_t varCount)
Returns whether for every variable .
Definition: Term.h:446
static void encodedDual(Exponent *res, const Exponent *dualOf, const Exponent *point, size_t varCount)
The parameter dualOf is interpreted to encode an irreducible ideal, and the dual of that reflected in...
Definition: Term.h:505
void displayNote(const string &msg)
Display msg to standard error in a way that indicates that this is something that the user should tak...
Definition: display.cpp:135
This file contains functions for printing strings to standard error.
void reportError(const string &errorMsg)
Definition: error.cpp:23
void codimension(const Ideal &ideal, mpz_t codim)
Compute the codimension of a monomial ideal.
Definition: frobby.cpp:441
void lcm(Word *res, const Word *resEnd, const Word *a, const Word *b)
#define ASSERT(X)
Definition: stdinc.h:86