1 #ifndef __ALGORITHM_BUCHBERGER_EXPLORER_CPP_ 2 #define __ALGORITHM_BUCHBERGER_EXPLORER_CPP_ 21 #include "hilbert_functions.hpp" 23 #include "algorithm_buchberger_basic.hpp" 25 #include "algorithm_buchberger_explorer.hpp" 27 #include "reduction_support.hpp" 47 if (not HPdiff.is_zero())
48 result = (HPdiff.numerator(HPdiff.degree()) < 0);
55 i <= h1->
degree() and i <= h2->
degree() and (*h1)[i] == (*h2)[i];
63 while (i < n and t[i] == u[i]) ++i;
64 if (i == n) result =
false;
72 if (i > h2->
degree()) result =
false;
73 else result = (*h1)[i] < (*h2)[i];
83 const list<Abstract_Polynomial *> &F,
86 WT_TYPE * strategy_weights,
87 const unsigned number_to_advance
89 unsigned number_of_spolys = 0;
90 list<Abstract_Polynomial *> G;
91 list<Critical_Pair_Basic *> P;
97 case StrategyFlags::NORMAL_STRATEGY:
break;
98 case StrategyFlags::SUGAR_STRATEGY:
101 case StrategyFlags::WSUGAR_STRATEGY:
110 bool verbose =
false;
111 bool very_verbose =
false;
127 for (; !P.empty() and i < number_to_advance; ++i) {
129 report_front_pair(Pp[i], strategy);
132 if ((s = Pp[i]->s_polynomial()) ==
nullptr)
140 while (i < number_to_advance) { Pp[i] =
nullptr; ++i; }
141 unsigned winning_index = number_to_advance + 1;
142 for (i = 0; i < number_to_advance and Pp[i] !=
nullptr; ++i) {
143 if (Pp[i] !=
nullptr) {
146 cout <<
"\treduced to zero\n";
148 HFN[i] = HSN[i]=
nullptr;
151 unsigned n = T.front().num_vars();
157 if (winning_index > number_to_advance)
160 if (newer_HF_smaller(
161 Pp[winning_index]->s_polynomial()->leading_monomial(),
163 Pp[i]->s_polynomial()->leading_monomial(),
172 cout <<
"sorting results\n";
173 for (i = 0; i < number_to_advance and Pp[i] !=
nullptr; ++i) {
174 delete HP[i];
delete HFN[i];
delete HSN[i];
175 if (i != winning_index or Pp[i]->s_polynomial()->is_zero()) {
178 if (not Pp[i]->s_polynomial()->is_zero())
182 cout <<
"\treduced to zero\n";
192 << ((Pp[i]->
is_generator()) ? 0 : Pp[i]->second()->leading_monomial())
194 if (Pp[i]->is_generator())
delete Pp[i]->first();
204 if (very_verbose) { cout <<
"\tadded "; r->println(); }
210 delete [] HP;
delete [] HFN;
delete [] HSN;
211 cout << number_of_spolys <<
" s-polynomials computed and reduced\n";
213 cout << G.
size() <<
" polynomials before interreduction\n";
216 cout << G.size() <<
" polynomials after interreduction\n";
218 list<Constant_Polynomial *> B;
The general class of a polynomial.
virtual bool is_zero() const =0
is this polynomial zero?
A Constant_Polynomial is a polynomial that should not change.
virtual void set_spoly(Mutable_Polynomial *p)
sets the s-polynomial to p, for explorer methods
void sort_pairs_by_strategy(list< T *> &P)
Applies the strategy to find the “smallest” critical pair.
Dense_Univariate_Rational_Polynomial * hilbert_polynomial(NVAR_TYPE n, unsigned int pole_order, const list< Monomial > T, Dense_Univariate_Integer_Polynomial *hn, Dense_Univariate_Integer_Polynomial *hn2)
computes the Hilbert polynomial for an ideal
void gm_update(list< Critical_Pair_Basic *> &P, list< Abstract_Polynomial *> &G, Abstract_Polynomial *r, StrategyFlags strategy)
Implementation of Gebauer-Moeller algorithm. Based on description in Becker and Weispfenning (1993)...
list< Constant_Polynomial * > buchberger_explorer(const list< Abstract_Polynomial *> &F, SPolyCreationFlags method, StrategyFlags strategy, WT_TYPE *strategy_weights, const unsigned number_to_advance)
Alternate implementation of Buchberger’s algorithm, for parallelization.
virtual Monomial & leading_monomial() const =0
leading monomial – call after sort_by_order()!
Dense_Univariate_Integer_Polynomial * hilbert_second_numerator(NVAR_TYPE n, Dense_Univariate_Integer_Polynomial *first, const WT_TYPE *grading)
computes the second Hilbert numerator (after reduction by )
Dense_Univariate_Integer_Polynomial * hilbert_numerator_bigatti(const list< Monomial > &T, const WT_TYPE *grading)
the Bigatti algorithm to compute the Hilbert numerator
quick-’n-dirty Dense_Univariate integer polynomial class
void reduce_over_basis(Mutable_Polynomial **sp, const T &G, int comm_id=0)
Reduce the polynomial r over the basis G.
virtual void at_generation_tasks()
hook called while first generating polynomial
StrategyFlags
flag indicating which strategy to use for computation
list< Abstract_Polynomial * > reduce_basis(list< Abstract_Polynomial *>G)
Remove redundant polynomials from G.
SPolyCreationFlags
flag indicating which structure to use for an s-polynomial
DEG_TYPE size
number of slots for coefficients
NVAR_TYPE num_vars() const
number of variables
Polynomials that need arithmetic typically descend from this class.
const Abstract_Polynomial * first() const
first polynomial in pair
polynomial-related data for a weighted sugar strategy
Implementation of monomials.
polynomial-related data for a sugar strategy
virtual Mutable_Polynomial * s_polynomial()
to use only if s-polynomial is already computed by another method
DEG_TYPE degree(NVAR_TYPE i) const
Degree of th variable.
void set_strategy(Poly_Strategy_Data *psd)
sets the polynomial’s strategy to psd
DEG_TYPE degree() const
the polynomial’s degree (exponent of largest nonzero term)
quick-’n-dirty Dense_Univariate rational polynomial class
Controls the creation of s-polynomials.
void report_critical_pairs(const list< T *>P, bool verbose=false)
A brief report on the number of critical pairs. If verbose is true, also lists them.
unsigned ideal_dimension(NVAR_TYPE n, const Dense_Univariate_Integer_Polynomial *h1, const Dense_Univariate_Integer_Polynomial *h2)
computes the dimension of the ideal by subtracting the Hilbert numerators
bool is_generator() const
whether this pair is from a generator
virtual Poly_Strategy_Data * strategy() const
strategy related information