1 #ifndef __ALGORITHM_BUCHBERGER_EXPLORER_CPP_ 2 #define __ALGORITHM_BUCHBERGER_EXPLORER_CPP_ 7 #include "hilbert_functions.hpp" 9 #include "algorithm_buchberger_basic.hpp" 10 #include "reduction_support.hpp" 12 #include "algorithm_buchberger_explorer.hpp" 16 #define XPLOR_GENERATOR -1 17 #define PROC_UNASSIGNED -1 32 pi = i;
qi = XPLOR_GENERATOR;
proc = PROC_UNASSIGNED;
36 int i,
int j,
unsigned strategy, vector<Abstract_Polynomial *>G
40 pi = i;
qi = j;
proc = PROC_UNASSIGNED;
45 vector<Abstract_Polynomial *>G
49 pi = i;
qi = G.size();
proc = PROC_UNASSIGNED;
60 return (static_cast<Pair_Sugar_Data *>(
key))->pair_sugar();
88 bool second_HF_smaller(
109 i <= h1->
degree() and i <= h2->
degree() and (*h1)[i] == (*h2)[i];
122 if (i > h2->
degree()) result =
false;
123 else result = (*h1)[i] < (*h2)[i];
156 i <= h1->
degree() and i <= h2->
degree() and (*h1)[i] == (*h2)[i];
164 while (i < n and t[i] == u[i]) ++i;
165 if (i == n) result =
false;
173 if (i > h2->
degree()) result =
false;
174 else result = (*h1)[i] < (*h2)[i];
195 void hilbert_chooser(
196 void * firstin,
void * secondin,
int * len, MPI_Datatype * dptr
198 int64_t *
first = (int64_t *)firstin;
199 int64_t *
second = (int64_t *)secondin;
200 bool undecided =
true;
201 bool keep_second =
true;
202 if (first[0] == -1 or second[0] == -1) {
209 int64_t n = first[1];
210 for (i = 0; undecided and i <= n; ++i)
211 undecided = (first[2*i + 2]*second[2*i+3] == second[2*i + 2]*first[2*i+3]);
214 keep_second = (first[2*i + 2]*second[2*i+3] >= second[2*i + 2]*first[2*i+3]);
217 for (i = 2*(n + 1) + 2; undecided and i < *len; ++i)
218 undecided = (first[i] == second[i]);
219 keep_second = undecided or first[--i] < second[i];
224 for (
unsigned i = 0; i < *len; ++i)
225 second[i] = first[i];
240 list<Critical_Pair_XPlor *> & P,
241 list<Critical_Pair_XPlor *> & Pass,
242 list<Critical_Pair_XPlor *> * Pcancel,
243 vector<Abstract_Polynomial *> & G,
248 list<Critical_Pair_XPlor *> C;
253 unsigned m = G.size();
254 for (
unsigned i = 0; i < G.size(); ++i)
258 list<Critical_Pair_XPlor *> D;
259 while (C.size() != 0) {
275 list<Critical_Pair_XPlor *> E;
276 while (D.size() != 0) {
290 list<Critical_Pair_XPlor *> Q;
291 while (P.size() != 0) {
307 list <Critical_Pair_XPlor *> R;
308 while (Pass.size() != 0) {
325 for (
unsigned i = 0; i < 2; ++i) {
326 if (Pcancel[i].size() > 0) {
327 cout <<
"Should delete pairs ";
329 cout <<
'(' <<
p->first_index() <<
',' <<
p->second_index() <<
" (" << i <<
") ) ";
354 const vector<Abstract_Polynomial *> &F,
357 WT_TYPE * strategy_weights,
361 double r_bcast_time = 0;
362 unsigned number_of_spolys = 0;
363 vector<Abstract_Polynomial *> G;
364 list<Critical_Pair_XPlor *> P;
365 list<Critical_Pair_XPlor *> Pass, * Pdel;
366 list<Critical_Pair_XPlor *> R;
367 list<Constant_Polynomial *> B;
371 MPI_Datatype pair_type, pair_basetypes[2];
372 MPI_Aint pair_offsets[2], extent, lb;
373 int pair_block_counts[2];
374 pair_offsets[0] = 0; pair_basetypes[0] = MPI_INT; pair_block_counts[0] = 2;
375 MPI_Type_get_extent(MPI_INT, &lb, &extent);
376 pair_offsets[1] = 2 * extent; pair_basetypes[1] = MPI_UNSIGNED_LONG_LONG;
377 pair_block_counts[1] = 1;
378 MPI_Type_create_struct(2, pair_block_counts, pair_offsets, pair_basetypes, &pair_type);
379 MPI_Type_commit(&pair_type);
381 MPI_Op hilbert_comparison_op;
382 MPI_Op_create(hilbert_chooser,
false, &hilbert_comparison_op);
384 vector<Abstract_Polynomial *>::const_iterator Fi = F.begin();
385 NVAR_TYPE n = (*Fi)->number_of_variables();
386 for (
unsigned i = 0; i < F.size(); ++i, ++Fi)
396 bool verbose =
false;
397 bool very_verbose =
false;
398 unsigned min_todo = 0;
400 Pdel =
new list<Critical_Pair_XPlor *>[comm_size];
409 MPI_Bcast(&min_todo, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD);
414 while (min_todo != 0) {
427 cout <<
"estimate " << min_todo <<
" steps left\n";
431 for (; !P.empty() and i < comm_size; ++i) {
439 cout << comm_id <<
" sending " 441 <<
" to " << i << endl;
442 MPI_Send(&p_new, 1, pair_type, i, 0, MPI_COMM_WORLD);
444 for (; P.empty() and i < comm_size; ++i) {
446 cout << comm_id <<
" sending " 448 <<
" to " << i << endl;
449 MPI_Send(&p_new, 1, pair_type, i, 0, MPI_COMM_WORLD);
459 cout << comm_id <<
" sending " 461 <<
" to " << comm_id << endl;
468 MPI_Recv(&p_in, 1, pair_type, 0, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
469 if (p_in.
first != XPLOR_GENERATOR) {
470 if (p_in.
second == XPLOR_GENERATOR)
490 if ((s =
p->s_polynomial()) ==
nullptr) {
491 s =
p->s_polynomial(method, strategy);
511 unsigned winning_index = R.size() + 1;
519 list<Critical_Pair_XPlor *>::iterator Rbest = R.end();
522 list<Critical_Pair_XPlor *>::iterator Ri = R.begin();
526 s = (*Ri)->s_polynomial();
528 cout <<
'\t' << comm_id <<
' ' << (*Ri)->first_index() <<
',' 529 << (*Ri)->second_index() <<
": reduced to zero\n";
531 if (Ri == R.begin()) {
535 list<Critical_Pair_XPlor *>::iterator Rdel = Ri;
542 unsigned n = T.front().num_vars();
551 if (second_HF_smaller(wHP, wHS, hp, hn)) {
581 int hp_deg = (wt ==
nullptr) ? -1 : hp->
degree();
583 MPI_Allreduce(&hp_deg, &recv_hp_deg, 1, MPI_INT, MPI_MAX, MPI_COMM_WORLD);
585 if (recv_hp_deg != -1) {
586 int hn_deg = (wHP ==
nullptr) ? -1 : wHS->
degree();
588 MPI_Allreduce(&hn_deg, &recv_hn_deg, 1, MPI_INT, MPI_MAX, MPI_COMM_WORLD);
589 int xfer_size = 2*(recv_hp_deg +1) + (recv_hn_deg + 1) + 2;
590 int64_t * send_data =
new int64_t[xfer_size];
594 send_data[0] = comm_id; send_data[1] = recv_hp_deg;
596 for ( ; hi < recv_hp_deg - wHP->
degree(); ++hi)
597 send_data[2*hi + 2] = send_data[2*hi + 3] = 0;
598 for ( ; hi <= recv_hp_deg; ++hi) {
599 send_data[2*hi + 2] = wHP->
numerator(recv_hp_deg - hi);
600 send_data[2*hi + 3] = wHP->
denominator(recv_hp_deg - hi);
602 for (hi = 0; hi <= wHS->
degree(); ++hi)
603 send_data[2 + 2*(recv_hp_deg + 1) + hi] = wHS->
coeff(hi);
605 hi = 2 + 2*(recv_hp_deg + 1) + wHS->
degree() + 1;
612 int64_t * recv_data =
new int64_t[xfer_size];
614 send_data, recv_data, xfer_size, MPI_INT64_T, hilbert_comparison_op,
617 winners_id = recv_data[0];
618 delete [] send_data;
delete [] recv_data;
624 if (winners_id >= 0) {
626 uint64_t r_bcast_size;
627 uint64_t r_bcast_sugar;
629 if (comm_id == winners_id) {
632 pair_data[0] = (*Rbest)->first_index();
633 pair_data[1] = (*Rbest)->second_index();
634 MPI_Send(pair_data, 2, MPI_INT, 0, 0, MPI_COMM_WORLD);
636 s = (*Rbest)->s_polynomial();
638 cout <<
"selected " << (*Rbest)->first_index()
653 MPI_Recv(pair_data, 2, MPI_INT, winners_id, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
655 list<Critical_Pair_XPlor *>::iterator Ri;
656 for (Ri = Pass.begin(); not found and Ri != Pass.end(); ) {
657 if ((*Ri)->first_index() == pair_data[0] and (*Ri)->second_index() == pair_data[1])
663 cout <<
"0 deleting " << pair_data[0] <<
',' << pair_data[1] << endl;
669 MPI_Bcast(&r_bcast_size, 1, MPI_UINT64_T, winners_id, MPI_COMM_WORLD);
670 if (comm_id != winners_id)
671 r_bcast =
new uint64_t [r_bcast_size];
673 MPI_Bcast(r_bcast, r_bcast_size, MPI_UINT64_T, winners_id, MPI_COMM_WORLD);
675 MPI_Bcast(&r_bcast_sugar, 1, MPI_UINT64_T, winners_id, MPI_COMM_WORLD);
677 if (comm_id != winners_id) {
678 r_bcast_size /= n + 1;
688 double start_time = MPI_Wtime();
691 very_verbose =
false;
692 if (very_verbose) { cout <<
"\tin full "; r->println(); }
693 very_verbose =
false;
696 for (
int i = 1; i < comm_size; ++i) {
697 int outgoing = Pdel[i].size();
698 MPI_Send(&outgoing, 1, MPI_INT, i, 0, MPI_COMM_WORLD);
699 while (Pdel[i].size() > 0) {
705 MPI_Send(&p_new, 1, pair_type, i, 0, MPI_COMM_WORLD);
710 while (Pdel[0].size() > 0) {
713 list<Critical_Pair_XPlor *>::iterator Ri = R.begin();
729 MPI_Recv(&incoming, 1, MPI_INT, 0, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
730 for (
int j = 0; j < incoming; ++j) {
731 MPI_Recv(&p_in, 1, pair_type, 0, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
732 list<Critical_Pair_XPlor *>::iterator Ri = R.begin();
735 ((*Ri)->first_index() != p_in.
first or
736 (*Ri)->second_index() != p_in.
second)
748 r_bcast_time += MPI_Wtime() - start_time;
763 if (comm_id == 0) min_todo = P.size();
764 unsigned my_todo = R.size();
766 MPI_Reduce(&my_todo, &all_todo, 1, MPI_UNSIGNED, MPI_MAX, 0, MPI_COMM_WORLD);
768 min_todo = min_todo < all_todo ? all_todo : min_todo;
769 MPI_Bcast(&min_todo, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD);
784 MPI_Type_free(&pair_type);
785 cout << comm_id <<
" reduced " << number_of_spolys << endl;
788 MPI_Reduce(&number_of_spolys, &total_spolys, 1, MPI_INT, MPI_SUM, 0, MPI_COMM_WORLD);
789 double max_bcast_time;
790 MPI_Reduce(&r_bcast_time, &max_bcast_time, 1, MPI_DOUBLE, MPI_MAX, 0, MPI_COMM_WORLD);
791 for (
unsigned i = 0; i < comm_size; ++i) {
793 cout << comm_id <<
" took " << r_bcast_time <<
" milliseconds on timed segment(s).\n";
794 MPI_Barrier(MPI_COMM_WORLD);
797 cout << total_spolys <<
" s-polynomials computed and reduced\n";
798 cout << max_bcast_time <<
" seconds spent on timed segment(s)\n";
800 cout << G.size() <<
" polynomials before interreduction\n";
802 list<Abstract_Polynomial *> G_final;
804 G_final.push_back(g);
806 cout << G_final.size() <<
" polynomials after interreduction\n";
Critical_Pair_XPlor(int i, Abstract_Polynomial *g, unsigned strategy, vector< Abstract_Polynomial *>G)
create critical pair (f,g) where f is at index i
The general class of a polynomial.
list< Constant_Polynomial * > buchberger_explorer(const vector< Abstract_Polynomial *> &F, SPolyCreationFlags method, StrategyFlags strategy, WT_TYPE *strategy_weights, const int comm_id, const int comm_size)
Alternate implementation of Buchberger’s algorithm, for parallelization.
virtual bool is_zero() const =0
is this polynomial zero?
MPZCOEF_TYPE numerator(DEG_TYPE k) const
returns the th numerator
int second
index of second polynomial in the pair
A Constant_Polynomial is a polynomial that should not change.
void sort_pairs_by_strategy(list< T *> &P)
Applies the strategy to find the “smallest” critical pair.
int first
index of first polynomial in the 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
bool is_zero() const
indicates whether the polynomial is zero
Critical_Pair_XPlor(int i, unsigned strategy, Abstract_Polynomial *f)
create critical pair (f,0) where f is at index i
const Abstract_Polynomial * second() const
second polynomial in pair
Mutable_Polynomial * s
S-polynomial.
void gm_update_explorer(list< Critical_Pair_XPlor *> &P, list< Critical_Pair_XPlor *> &Pass, list< Critical_Pair_XPlor *> *Pcancel, vector< Abstract_Polynomial *> &G, Abstract_Polynomial *r, unsigned strategy)
Implementation of Gebauer-Moeller algorithm, with XPLOR critical pairs. Based on description in Becke...
Critical_Pair_XPlor(int i, int j, unsigned strategy, vector< Abstract_Polynomial *>G)
create critical pair (f,g) where f, g are at indices i, j
int first_index()
returns index of first polynomial in pair
used to pass inforation on a critical pair from one polynomial to another
virtual Monomial & leading_monomial() const =0
leading monomial – call after sort_by_order()!
const Monomial & lcm() const
lcm of leading monomials of polynomials
int get_processor()
query whether this pair is assigned to a processor, and which (nonnegative value indicates assignment...
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 )
void set_processor(int i)
record that this pair is assigned to processor i
int second_index()
returns index of second polynomial in pair
COEF_TYPE coeff(DEG_TYPE k) const
the th coefficient
uint64_t * serialized(uint64_t &size)
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
DEG_TYPE poly_sugar() const
the sugar itself
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
list< Abstract_Polynomial * > reduce_basis(list< Abstract_Polynomial *>G)
Remove redundant polynomials from G.
DEG_TYPE size
number of slots for coefficients
int qi
second polynomial in the critical pair
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.
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
virtual Monomial & leading_monomial() const override
Implementation of monomials.
interface to a monomial ordering
polynomial-related data for a sugar strategy
Encapsulates information about a polynomial ring for easy access: ground field, number of indetermina...
int pi
first polynomial in the critical pair
Critical_Pair_XPlor(int i, StrategyFlags strategy, Abstract_Polynomial *f)
create critical pair (f,0) where f is at index i
bool is_coprime(const Monomial &other) const
true iff this has no common factor with other.
DEG_TYPE degree(NVAR_TYPE i) const
Degree of th variable.
Pair_Strategy_Data * key
strategy used to sort critical pairs
void set_strategy(Poly_Strategy_Data *psd)
sets the polynomial’s strategy to psd
DEG_TYPE sugar()
returns sugar of this pair; use ONLY if with sugar strategy
void force_sugar(DEG_TYPE new_sugar)
for those times when a different sugar is appropriate
DEG_TYPE degree() const
the polynomial’s degree (exponent of largest nonzero term)
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
int proc
processor to which this pair has been assigned
MPZCOEF_TYPE denominator(DEG_TYPE k) const
returns the th denominator
DEG_TYPE degree() const
returns the polynomial’s degree
Controls the creation of s-polynomials.
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
DEG_TYPE sugar
pair’s sugar
Abstract_Polynomial * p
first polynomial in the critical pair
contains information on critical pairs by their index in the basis, in addition to the usual informat...
virtual Poly_Strategy_Data * strategy() const
strategy related information