Gröbner basis project
Codebase for research into Gröbner basis computation
algorithm_buchberger_dynamic.cpp
1 #ifndef __ALGORITHM_BUCHBERGER_DYNAMIC_CPP_
2 #define __ALGORITHM_BUCHBERGER_DYNAMIC_CPP_
3 
4 /*****************************************************************************\
5 * This file is part of DynGB. *
6 * *
7 * DynGB is free software: you can redistribute it and/or modify *
8 * it under the terms of the GNU General Public License as published by *
9 * the Free Software Foundation, either version 2 of the License, or *
10 * (at your option) any later version. *
11 * *
12 * Foobar is distributed in the hope that it will be useful, *
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
15 * GNU General Public License for more details. *
16 * *
17 * You should have received a copy of the GNU General Public License *
18 * along with DynGB. If not, see <http://www.gnu.org/licenses/>. *
19 \*****************************************************************************/
20 
21 #include "system_constants.hpp"
22 
23 #include "algorithm_buchberger_basic.hpp"
24 #include "algorithm_buchberger_dynamic.hpp"
25 
26 #include "skeleton.hpp"
27 #include "dynamic_engine.hpp"
28 #include "particular_orderings.hpp"
29 #include "reduction_support.hpp"
30 
31 // caching the weights seems to run about 15% slower, oddly enough
32 //#define ORDERING_TYPE CachedWGrevlex_Ordering
33 #define ORDERING_TYPE WGrevlex
34 
36  Mutable_Polynomial **sp,
37  const list<Abstract_Polynomial *>G
38 ) {
39  Abstract_Polynomial * g; // used to loop through G
40  Mutable_Polynomial * s = *sp; // s-poly
41  Mutable_Polynomial * r = s->zero_polynomial(); // remainder / residue
42  bool verbose = false;
43  bool very_verbose = false;
44  // continue reducing until s is zero
45  while (!s->is_zero()) {
46  if (verbose) cout << "reducing " << s->leading_monomial() << "? ";
47  if (very_verbose) s->println();
48  if ((g = find_reducer<list<Abstract_Polynomial *>>(s, G))) {
49  if (verbose) cout << "\tyes! by " << g->leading_monomial() << endl;
50  if (very_verbose) { cout << "\tyes! by "; g->println(); }
51  // DYNAMIC
52  if (g->monomial_ordering() != s->monomial_ordering())
54  Monomial t(s->leading_monomial());
55  t /= g->leading_monomial();
58  a /= b;
59  // s = s - atg
60  if (s->strategy() != nullptr)
61  s->strategy()->pre_reduction_tasks(t, *g);
62  s->add_polynomial_multiple(a, t, *g, true);
63  if (very_verbose) {
64  cout << "\tresults in ";
65  if (s->is_zero())
66  cout << 0 << endl;
67  else
68  s->println();
69  }
70  if (verbose and not s->is_zero())
71  cout << "\tresults in " << s->leading_monomial() << endl;
72  } else {
73  if (verbose or very_verbose) cout << "no!\n";
74  // no reducer; move leading term to residue and continue
76  // t should already be smaller than what's already in r, so no danger
78  delete t;
79  }
80  }
81  // move the strategy from s to r
82  r->set_strategy(s->strategy());
83  s->set_strategy(nullptr);
84  delete s;
85  *sp = r;
86 }
87 
100  list<Critical_Pair_Dynamic *> & P,
101  list<Abstract_Polynomial *> & G,
103  StrategyFlags strategy,
104  ORDERING_TYPE * ordering
105 ) {
106  //cout << "----------------------\n";
107  list<Critical_Pair_Dynamic *> C;
108  // critical pairs with new polynomial
109  for (Abstract_Polynomial * g : G)
110  C.push_back(new Critical_Pair_Dynamic(g, r, strategy, ordering));
111  // apply Buchberger's lcm criterion to new pairs
112  list<Critical_Pair_Dynamic *> D;
113  while (C.size() != 0) {
114  Critical_Pair_Dynamic * p = C.front();
115  C.pop_front();
116  if ((p->first()->leading_monomial().is_coprime(
117  p->second()->leading_monomial()))
118  or (no_triplet(p, C) and no_triplet(p, D))
119  )
120  D.push_back(p);
121  else {
122  //cout << "triplet prunes " << *p << endl;
123  delete p;
124  }
125  }
126  // apply Buchberger's gcd criterion
127  list<Critical_Pair_Dynamic *> E;
128  while (D.size() != 0) {
129  Critical_Pair_Dynamic * p = D.front();
130  D.pop_front();
131  if (!(p->first()->leading_monomial().is_coprime(
132  p->second()->leading_monomial())))
133  E.push_back(p);
134  else {
135  //cout << "gcd prunes " << *p << endl;
136  delete p;
137  }
138  }
139  // apply Buchberger's lcm criterion to old pairs
140  list<Critical_Pair_Dynamic *> Q;
141  while (P.size() != 0) {
142  Critical_Pair_Dynamic * p = P.front();
143  P.pop_front();
144  if (!(r->leading_monomial() | p->lcm())
145  or lcm_alike(p->first()->leading_monomial(), r->leading_monomial(), p)
147  )
148  Q.push_back(p);
149  else {
150  //cout << "triplet prunes " << *p << endl;
151  delete p;
152  }
153  }
154  P = Q;
155  // add new pairs to old pairs
156  for (Critical_Pair_Dynamic * e : E)
157  P.push_back(e);
158  /*cout << "All pairs:\n";
159  for (auto pi = P.begin(); pi != P.end(); ++pi)
160  cout << '\t' << **pi << endl;
161  cout << "----------------------\n";*/
162  // add new poly to basis
163  G.push_back(r);
164 }
165 
166 bool none_others_divide(const list<Monomial> & U, const Monomial & t) {
167  bool result = true;
168  for (const Monomial & u : U) {
169  if (not u.is_like(t) and u | t) {
170  result = false;
171  break;
172  }
173  }
174  return result;
175 }
176 
177 list<list<Monomial>> recursive_select(list<list<Monomial>> & Ms) {
178  list<list<Monomial>> result;
179  if (Ms.size() == 0) {
180  list<Monomial> M;
181  result.push_back(M);
182  } else {
183  auto M = Ms.front();
184  Ms.pop_front();
185  auto lower = recursive_select(Ms);
186  for (auto t: M) {
187  if (none_others_divide(M, t)) {
188  for (auto S : lower) {
189  auto tmp(S);
190  tmp.push_front(t);
191  result.push_back(tmp);
192  }
193  }
194  }
195  }
196  return result;
197 }
198 
199 void initial_analysis(
200  const list<Abstract_Polynomial *> & F,
201  Monomial_Ordering * & mord,
202  LP_Solver * solver
203 ) {
204  NVAR_TYPE n = (*F.begin())->leading_monomial().num_vars();
205  WT_TYPE weights[n];
206  // arrange possible monomial chains
207  list<list<Monomial>> mons;
208  list<list<Monomial>> saved_mons;
209  for (auto f : F) {
210  set<Monomial> T;
211  set<Monomial> U;
212  set<Monomial> V;
213  Polynomial_Iterator * fi = f->new_iterator();
214  while (not fi->fellOff()) {
215  U.insert((fi->currMonomial()));
216  fi->moveRight();
217  }
218  delete fi;
219  compatiblePP(f->leading_monomial(), U, solver->get_rays(), T, V, solver);
220  list<Monomial> W;
221  for (auto t : T) W.push_back(t);
222  mons.push_back(W);
223  saved_mons.push_back(W);
224  }
225  list<list<Monomial>> mon_chains = recursive_select(mons);
226  // now check each chain for solvability & desirability
227  Dense_Univariate_Integer_Polynomial * win_hn = nullptr;
228  Dense_Univariate_Rational_Polynomial * win_hp = nullptr;
229  vector<constraint> all_constraints;
230  for (auto M : mon_chains) {
231  list<Monomial> T;
232  list<Abstract_Polynomial *> G;
233  list<Critical_Pair_Dynamic *> P;
236  = hilbert_second_numerator(n, T_hn);
238  = hilbert_polynomial(n, ideal_dimension(n, T_hn, T_hn2), M, T_hn, T_hn2);
239  delete T_hn2;
240  if (not (win_hn == nullptr or hilbertCmp(*win_hn, *win_hp, *T_hn, *T_hp) > 0))
241  {
242  delete T_hn;
243  delete T_hp;
244  } else {
245  // temporarily solve using GLPK
246  GLPK_Solver * tmp_lp = new GLPK_Solver(n);
247  set<constraint> tmp_cnstr;
248  bool consistent = true;
249  CONSTR_TYPE w[n];
250  auto Mi = M.begin();
251  for (auto N : saved_mons) {
252  auto t = *Mi;
253  auto fi = N.begin();
254  while (consistent and fi != N.end()) {
255  if (not t.is_like(*fi)) {
256  for (NVAR_TYPE i = 0; i < n; ++i)
257  w[i] = t.degree(i) - fi->degree(i);
258  constraint c(n, w);
259  consistent = tmp_lp->solve(c);
260  if (consistent) tmp_cnstr.emplace(c);
261  }
262  ++fi;
263  }
264  if (not consistent)
265  break;
266  ++Mi;
267  }
268  delete tmp_lp;
269  if (consistent) {
270  all_constraints.clear();
271  for (auto c : tmp_cnstr)
272  all_constraints.push_back(c);
273  if (win_hn != nullptr) {
274  delete win_hn;
275  delete win_hp;
276  }
277  win_hn = T_hn;
278  win_hp = T_hp;
279  }
280  }
281  }
282  solver->solve(all_constraints);
283  cout << "prefer ordering with " << solver->get_rays().size() << " rays: ";
284  ::ray interior_ray = ray_sum(solver->get_rays());
285  cout << interior_ray << endl;
286  auto wt_ptr = interior_ray.weights();
287  for (NVAR_TYPE i = 0; i < n; ++i) weights[i] = wt_ptr[i];
288  if (mord != nullptr)
289  delete mord;
290  mord = new ORDERING_TYPE(n, weights);
291 }
292 
293 list<Constant_Polynomial *> buchberger_dynamic(
294  const list<Abstract_Polynomial *> &F,
295  SPolyCreationFlags method,
296  StrategyFlags strategy,
297  WT_TYPE * strategy_weights,
298  DynamicHeuristic heuristic,
299  DynamicSolver solver_type,
300  bool analyze_inputs
301 ) {
302  unsigned number_of_spolys = 0;
303  list<Abstract_Polynomial *> G;
304  list<Critical_Pair_Dynamic *> P;
305  // set up list of LPPs for selecting monomials
306  list<Monomial> currentLPPs;
307  // create initial ordering
308  NVAR_TYPE n = (*(F.begin()))->base_ring().number_of_variables();
309  WT_TYPE * weights = new WT_TYPE [n];
310  for (NVAR_TYPE i = 0; i < n; ++i) weights[i] = 1;
311  ORDERING_TYPE * curr_ord = new ORDERING_TYPE(n, weights);
312  list<ORDERING_TYPE *> all_orderings_used; // so we can free them at the end
313  all_orderings_used.push_front(curr_ord);
314  const WT_TYPE * wt_ptr;
315  // create skeleton
316  LP_Solver * solver;
317  switch (solver_type) {
318  case SKELETON_SOLVER: solver = new skeleton(n); break;
319  case GLPK_SOLVER: solver = new GLPK_Solver(n); break;
320  case PPL_SOLVER: solver = new PPL_Solver(n); break;
321  default: solver = new skeleton(n); break;
322  }
323  //G.push_back(*(F.begin()));
324  // set up initial ordering
325  Abstract_Polynomial * g = new Constant_Polynomial(**(F.begin()));
326  cout << "Working with " << g->leading_monomial() << endl;
327  switch (strategy) {
328  case StrategyFlags::SUGAR_STRATEGY:
329  g->set_strategy(new Poly_Sugar_Data(g));
330  break;
331  case StrategyFlags::WSUGAR_STRATEGY:
332  g->set_strategy(new Poly_WSugar_Data(g, strategy_weights));
333  break;
334  default: break; // includes NORMAL_STRATEGY
335  }
336  bool new_world_order = false;
337  Dense_Univariate_Integer_Polynomial * hNum = nullptr;
338  SelectMonomial(g, currentLPPs, &hNum, G, P, solver, new_world_order, heuristic);
339  G.push_back(g);
340  // check ordering & convert if necessary
341  if (new_world_order) {
342  weights = new WT_TYPE [n];
343  cout << "new ordering from " << solver->get_rays().size() << " rays: ";
344  ::ray interior_ray = ray_sum(solver->get_rays());
345  cout << interior_ray << endl;
346  wt_ptr = interior_ray.weights();
347  for (NVAR_TYPE i = 0; i < n; ++i) weights[i] = wt_ptr[i];
348  curr_ord = new ORDERING_TYPE(n, weights);
349  all_orderings_used.push_front(curr_ord);
350  g->set_monomial_ordering(curr_ord);
351  }
352  /*cout << "Concluded with " << g->leading_monomial() << endl;
353  for (auto r: solver->get_rays()) cout << '\t' << r << endl;
354  cout << endl;*/
355  // set up critical_pairs
356  for (Abstract_Polynomial * f : F)
357  if (f != *(F.begin())) {
358  Critical_Pair_Dynamic * p = new Critical_Pair_Dynamic(f, strategy, curr_ord);
359  P.push_back(p);
360  }
361  // main loop
362  bool verbose = false;
363  bool very_verbose = false;
364  while (!P.empty()) {
365  //cout << "current basis: ";
366  //for (auto g : G) cout << g->leading_monomial() << ',';
367  //cout << endl;
368  cout << "current ordering: " << *curr_ord << endl;
369  sort_pairs_by_strategy<Critical_Pair_Dynamic>(P);
371  Critical_Pair_Dynamic * p = P.front();
372  P.pop_front();
373  report_front_pair(p, strategy);
374  cout << "\tweighted degree: " << p->lcm().weighted_degree(weights) << endl;
375  /*cout << "\tothers:\n\t\t";
376  for (Critical_Pair_Dynamic * pp : P) cout << pp->lcm() << ',' << static_cast<const Pair_Sugar_Data *>(pp->pair_key())->pair_sugar() << ',' << pp->lcm().weighted_degree(weights) << "; ";
377  cout << endl;*/
378  // make s-poly
379  // DYNAMIC
380  Mutable_Polynomial * s = p->s_polynomial(method, strategy);
381  ++number_of_spolys;
382  delete p;
383  // cout << "Reducing s-poly "; s->println();
384  if (!s->is_zero())
386  if (s->is_zero()) {
387  cout << "\treduced to zero\n";
388  delete s;
389  } else {
391  r = new Constant_Polynomial(*s);
392  // move strategy from s to r
393  r->set_strategy(s->strategy());
394  s->set_strategy(nullptr);
395  delete s;
396  cout << "\tadded " << r->leading_monomial() << endl;
397  very_verbose = false;
398  if (very_verbose) { cout << "\tadded "; r->println(); }
399  very_verbose = false;
400  // Currently performing "late conversion" of polynomials to new orderings
401  // so new_world_order should be useless for the time being
402  bool new_world_order = false;
403  SelectMonomial(r, currentLPPs, &hNum, G, P, solver, new_world_order, heuristic);
404  if (new_world_order) {
405  cout << "new ordering from " << solver->get_rays().size() << " rays: ";
406  ::ray interior_ray = ray_sum(solver->get_rays());
407  cout << interior_ray << endl;
408  //for (auto r: solver->get_rays()) cout << '\t' << r << endl;
409  //cout << endl;
410  wt_ptr = interior_ray.weights();
411  weights = new WT_TYPE [n];
412  for (NVAR_TYPE i = 0; i < n; ++i) weights[i] = wt_ptr[i];
413  curr_ord = new ORDERING_TYPE(n, weights);
414  //cout << skel << endl;
415  all_orderings_used.push_front(curr_ord);
416  r->set_monomial_ordering(curr_ord);
417  cout << "ended with leading monomial " << r->leading_monomial() << endl;
418  //cout << "changing " << P.size() << " pairs:\n\t";
419  for (Critical_Pair_Dynamic * p : P) {
420  //cout << p->lcm() << " ";
421  p->change_ordering(curr_ord);
422  }
423  cout << endl;
424  }
425  gm_update_dynamic(P, G, r, strategy, curr_ord);
426  }
427  }
428  cout << number_of_spolys << " s-polynomials computed and reduced\n";
429  // cleanup
430  cout << G.size() << " polynomials before interreduction\n";
431  //check_correctness(G, strategy);
432  G = reduce_basis(G);
433  cout << G.size() << " polynomials after interreduction\n";
434  //set<Constant_Polynomial *, smaller_lm> B;
435  // sort all polynomials to new ordering
436  list<Constant_Polynomial *> B;
437  unsigned long num_mons = 0;
438  unsigned long max_mons = 0;
439  for (Abstract_Polynomial * g : G)
440  {
441  unsigned long glen = g->length();
442  num_mons += glen;
443  if (glen > max_mons) max_mons = glen;
445  b->set_monomial_ordering(curr_ord);
446  B.push_back(b);
447  delete g;
448  }
449  cout << "tot # monomials: " << num_mons << endl;
450  cout << "max # monomials: " << max_mons << endl;
451  cout << "avg # monomials: " << num_mons / B.size() << endl;
452  // eliminate old orderings
453  all_orderings_used.pop_front();
454  // first one should be curr_ord; we do not want to delete that!
455  while (all_orderings_used.size() != 0) {
456  ORDERING_TYPE * bye_bye_ordering = all_orderings_used.front();
457  all_orderings_used.pop_front();
458  delete [] bye_bye_ordering->order_weights();
459  delete bye_bye_ordering;
460  }
461  if (hNum != nullptr) delete hNum;
462  delete solver;
463  return B;
464 }
465 
466 #endif
The general class of a polynomial.
Definition: polynomial.hpp:101
ray ray_sum(const set< ray > &rs)
Add all the rays in a set.
Definition: lp_solver.cpp:281
virtual bool is_zero() const =0
is this polynomial zero?
void SelectMonomial(Abstract_Polynomial *r, list< Monomial > &CurrentLPPs, Dense_Univariate_Integer_Polynomial **current_hilbert_numerator, const list< Abstract_Polynomial *> &CurrentPolys, const list< Critical_Pair_Dynamic *> &crit_pairs, LP_Solver *currSkel, bool &ordering_changed, DynamicHeuristic method)
Selects a leading power product for a polynomial.
const Monomial_Ordering * monomial_ordering() const
reports leading monomial’s monomial ordering
Definition: polynomial.hpp:130
DEG_TYPE weighted_degree(const WT_TYPE *weights, NVAR_TYPE m=0) const
Definition: monomial.cpp:199
A Constant_Polynomial is a polynomial that should not change.
void change_ordering(Weighted_Ordering *new_order)
change the ordering associated with this pair
approximate skeleton of a polyhedral cone, using PPL linear solver
Definition: ppl_solver.hpp:38
void reduce_over_basis_dynamic(Mutable_Polynomial **sp, const list< Abstract_Polynomial *>G)
Reduce the polynomial r over the basis G.
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
virtual void moveRight()=0
Moves right in the polynomial, to the next smaller monomial.
const Abstract_Polynomial * second() const
second polynomial in pair
virtual Mutable_Polynomial * s_polynomial(SPolyCreationFlags method, StrategyFlags strategy)
creates the s-polynomial for first() and second()
virtual const set< ray > & get_rays()
Returns the rays that define the skeleton.
Definition: lp_solver.cpp:378
virtual Monomial & leading_monomial() const =0
leading monomial – call after sort_by_order()!
virtual Prime_Field_Element leading_coefficient() const =0
leading coefficient – call after sort_by_order()!
const Monomial & lcm() const
lcm of leading monomials of polynomials
virtual Mutable_Polynomial * zero_polynomial() const =0
zero polynomial of this type
Controls the creation of s-polynomials, specialized for the dynamic algorithm.
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 )
approximate skeleton of a polyhedral cone, using GLPK linear solver
Definition: glpk_solver.hpp:37
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
StrategyFlags
flag indicating which strategy to use for computation
Definition: strategies.hpp:34
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
a constraint
Definition: lp_solver.hpp:53
skeleton of a polyhedral cone, with methods allowing definition and refinement
Definition: skeleton.hpp:178
virtual void add_last(const Prime_Field_Element &, const Monomial &)=0
Attach a new monomial to the tail – check that it belongs at tail!
void compatiblePP(Monomial currentLPP, const set< Monomial > &allPPs, const set<::ray > &bndrys, set< Monomial > &result, set< Monomial > &boundary_mons, LP_Solver *skel)
DynamicHeuristic
bool lcm_alike(const Monomial &t, const Monomial &u, const Critical_Pair_Basic *p)
Checks if the lcm of t and u is like the lcm stored in p.
exact or approximate polyhedral cone solution, with methods allowing definition and refinement ...
Definition: lp_solver.hpp:504
virtual bool fellOff() const =0
Polynomials that need arithmetic typically descend from this class.
Definition: polynomial.hpp:305
const Abstract_Polynomial * first() const
first polynomial in pair
polynomial-related data for a weighted sugar strategy
Implementation of monomials.
Definition: monomial.hpp:69
Element of a field of prime characteristic.
Definition: fields.hpp:137
interface to a monomial ordering
virtual Mutable_Polynomial * detach_head()=0
Remove and return the head.
virtual void set_monomial_ordering(const Monomial_Ordering *order, bool sort_anew=true)=0
set the monomial ordering and sort the polynomials (optionally, but by default)
polynomial-related data for a sugar strategy
bool is_coprime(const Monomial &other) const
true iff this has no common factor with other.
Definition: monomial.cpp:243
int hilbertCmp(const Dense_Univariate_Integer_Polynomial &hn1, const Dense_Univariate_Rational_Polynomial &hp1, const Dense_Univariate_Integer_Polynomial &hn2, const Dense_Univariate_Rational_Polynomial &hp2)
DEG_TYPE degree(NVAR_TYPE i) const
Degree of th variable.
Definition: monomial.cpp:183
virtual void set_monomial_ordering(const Monomial_Ordering *order, bool sort_anew=true) override
sets the ordering of monomials in this polynomial
virtual void add_polynomial_multiple(const Prime_Field_Element &, const Monomial &, const Abstract_Polynomial &, bool subtract=false)=0
add monomial multiple of other
DynamicSolver
used by buchberger_dynamic() to decide which solver to use
void set_strategy(Poly_Strategy_Data *psd)
sets the polynomial’s strategy to psd
Definition: polynomial.cpp:36
list< Constant_Polynomial * > buchberger_dynamic(const list< Abstract_Polynomial *> &F, SPolyCreationFlags method, StrategyFlags strategy, WT_TYPE *strategy_weights, DynamicHeuristic heuristic, DynamicSolver solver_type, bool analyze_inputs)
implementation of the dynamic Buchberger algorithm
Used to iterate through a polynomial.
Definition: polynomial.hpp:224
void gm_update_dynamic(list< Critical_Pair_Dynamic *> &P, list< Abstract_Polynomial *> &G, Abstract_Polynomial *r, StrategyFlags strategy, ORDERING_TYPE *ordering)
implementation of Gebauer-Möller algorithm, adjusted for dynamic computation
virtual const Monomial & currMonomial() const =0
Reports the monomial at the current position.
bool no_triplet(const T *p, const list< T *>C)
Checks whether p is in danger of forming a Buchberger triple with some pair listed in C...
quick-’n-dirty Dense_Univariate rational polynomial class
bool is_like(const Monomial &other) const
Have same variables, same powers? Synonymous with operator==().
Definition: monomial.cpp:247
virtual bool solve(constraint &)=0
Adds the indicated constraint (singular!) and re-computes the solution.
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
const RAYENT_TYPE * weights() const
Returns the weights.
Definition: lp_solver.hpp:265
virtual unsigned length() const =0
number of monomials
virtual bool solve(constraint &)
Adds the indicated constraint (singular!) and re-computes the solution.
Abstract_Polynomial * find_reducer(Abstract_Polynomial *r, const T &G)
Find a polynomial in the basis G that can reduce r.
a ray defined by nonnegative coordinates
Definition: lp_solver.hpp:190
virtual Poly_Strategy_Data * strategy() const
strategy related information
Definition: polynomial.hpp:144