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];
194 void hilbert_chooser(
195 void * firstin,
void * secondin,
int * len, MPI_Datatype * dptr
197 int64_t *
first = (int64_t *)firstin;
198 int64_t *
second = (int64_t *)secondin;
199 bool undecided =
true;
200 bool keep_second =
true;
202 int64_t n = first[1];
203 for (i = 2; undecided and i < n + 2; ++i)
204 undecided = keep_second = (first[i] < second[i]);
206 for ( ; undecided and i < *len; ++i)
207 undecided = keep_second = (first[i] == second[i]);
210 for (
unsigned i = 0; i < *len; ++i)
211 second[i] = first[i];
226 list<Critical_Pair_XPlor *> & P,
227 list<Critical_Pair_XPlor *> & Pass,
228 list<Critical_Pair_XPlor *> * Pcancel,
229 vector<Abstract_Polynomial *> & G,
234 list<Critical_Pair_XPlor *> C;
239 unsigned m = G.size();
240 for (
unsigned i = 0; i < G.size(); ++i)
244 list<Critical_Pair_XPlor *> D;
245 while (C.size() != 0) {
261 list<Critical_Pair_XPlor *> E;
262 while (D.size() != 0) {
276 list<Critical_Pair_XPlor *> Q;
277 while (P.size() != 0) {
293 list <Critical_Pair_XPlor *> R;
294 while (Pass.size() != 0) {
332 const vector<Abstract_Polynomial *> &F,
335 WT_TYPE * strategy_weights,
339 double r_bcast_time = 0;
340 unsigned number_of_spolys = 0;
341 vector<Abstract_Polynomial *> G;
342 list<Critical_Pair_XPlor *> P;
343 list<Critical_Pair_XPlor *> Pass, * Pdel;
344 list<Critical_Pair_XPlor *> R;
345 list<Constant_Polynomial *> B;
349 MPI_Datatype pair_type, pair_basetypes[2];
350 MPI_Aint pair_offsets[2], extent, lb;
351 int pair_block_counts[2];
352 pair_offsets[0] = 0; pair_basetypes[0] = MPI_INT; pair_block_counts[0] = 2;
353 MPI_Type_get_extent(MPI_INT, &lb, &extent);
354 pair_offsets[1] = 2 * extent; pair_basetypes[1] = MPI_UNSIGNED_LONG_LONG;
355 pair_block_counts[1] = 1;
356 MPI_Type_create_struct(2, pair_block_counts, pair_offsets, pair_basetypes, &pair_type);
357 MPI_Type_commit(&pair_type);
359 MPI_Op hilbert_comparison_op;
360 MPI_Op_create(hilbert_chooser,
false, &hilbert_comparison_op);
362 vector<Abstract_Polynomial *>::const_iterator Fi = F.begin();
363 NVAR_TYPE n = (*Fi)->number_of_variables();
364 for (
unsigned i = 0; i < F.size(); ++i, ++Fi)
374 bool verbose =
false;
375 bool very_verbose =
false;
376 unsigned min_todo = 0;
378 Pdel =
new list<Critical_Pair_XPlor *>[comm_size];
387 MPI_Bcast(&min_todo, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD);
392 while (min_todo != 0) {
408 for (; !P.empty() and i < comm_size; ++i) {
419 MPI_Send(&p_new, 1, pair_type, i, 0, MPI_COMM_WORLD);
421 for (; P.empty() and i < comm_size; ++i) {
426 MPI_Send(&p_new, 1, pair_type, i, 0, MPI_COMM_WORLD);
445 MPI_Recv(&p_in, 1, pair_type, 0, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
446 if (p_in.
first != XPLOR_GENERATOR) {
447 if (p_in.
second == XPLOR_GENERATOR)
456 if ((s =
p->s_polynomial()) ==
nullptr) {
457 s =
p->s_polynomial(method, strategy);
474 double start_time = MPI_Wtime();
475 unsigned winning_index = R.size() + 1;
483 list<Critical_Pair_XPlor *>::iterator Rbest = R.end();
486 list<Critical_Pair_XPlor *>::iterator Ri = R.begin();
490 s = (*Ri)->s_polynomial();
495 if (Ri == R.begin()) {
499 list<Critical_Pair_XPlor *>::iterator Rdel = Ri;
506 unsigned n = T.front().num_vars();
515 if (second_HF_smaller(wHP, wHS, hp, hn)) {
535 int64_t size_data = -1;
537 if (wHP !=
nullptr) {
538 size_data = wHP->
degree();
539 MPI_Send(&size_data, 1, MPI_INT64_T, 0, 0, MPI_COMM_WORLD);
540 MPI_Send(wHP->
numerator_array(), size_data + 1, MPI_INT64_T, 0, 0, MPI_COMM_WORLD);
541 MPI_Send(wHP->
denominator_array(), size_data + 1, MPI_UINT64_T, 0, 0, MPI_COMM_WORLD);
543 size_data = wHS->
degree();
544 MPI_Send(&size_data, 1, MPI_INT64_T, 0, 0, MPI_COMM_WORLD);
545 MPI_Send(wHS->
coefficient_array(), size_data + 1, MPI_INT64_T, 0, 0, MPI_COMM_WORLD);
549 MPI_Send(&size_data, 1, MPI_INT64_T, 0, 0, MPI_COMM_WORLD);
564 if (wHP !=
nullptr) winners_id = 0;
565 for (
int i = 1; i < comm_size; ++i) {
566 MPI_Recv(&size_data, 1, MPI_INT64_T, i, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
567 if (size_data >= 0) {
568 int64_t * nums =
new int64_t[size_data + 1];
569 uint64_t * denoms =
new uint64_t[size_data + 1];
570 MPI_Recv(nums, size_data + 1, MPI_INT64_T, i, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
571 MPI_Recv(denoms, size_data + 1, MPI_UINT64_T, i, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
579 delete [] nums;
delete [] denoms;
580 MPI_Recv(&size_data, 1, MPI_UINT64_T, i, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
581 int64_t * coefs =
new int64_t[size_data + 1];
582 MPI_Recv(coefs, size_data + 1, MPI_INT64_T, i, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
587 if (wHP ==
nullptr) {
588 winners_id = i; wHP = hp_in; wHS = hn_in;
590 else if (second_HF_smaller(wHP, wHS, hp_in, hn_in)) {
591 if (wHP !=
nullptr)
delete wHP;
592 if (wHS !=
nullptr)
delete wHS;
593 winners_id = i; wHP = hp_in; wHS = hn_in;
595 delete hp_in;
delete hn_in;
602 if (wHP ==
nullptr) winners_id = -1;
609 MPI_Bcast(&winners_id, 1, MPI_INT, 0, MPI_COMM_WORLD);
610 if (winners_id >= 0) {
612 uint64_t r_bcast_size;
613 uint64_t r_bcast_sugar;
615 if (comm_id == winners_id) {
617 s = (*Rbest)->s_polynomial();
633 MPI_Bcast(&r_bcast_size, 1, MPI_UINT64_T, winners_id, MPI_COMM_WORLD);
634 if (comm_id != winners_id)
635 r_bcast =
new uint64_t [r_bcast_size];
637 MPI_Bcast(r_bcast, r_bcast_size, MPI_UINT64_T, winners_id, MPI_COMM_WORLD);
639 MPI_Bcast(&r_bcast_sugar, 1, MPI_UINT64_T, winners_id, MPI_COMM_WORLD);
641 if (comm_id != winners_id) {
642 r_bcast_size /= n + 1;
648 r_bcast_time += MPI_Wtime() - start_time;
652 very_verbose =
false;
653 if (very_verbose) { cout <<
"\tin full "; r->println(); }
654 very_verbose =
false;
657 for (
int i = 1; i < comm_size; ++i) {
658 int outgoing = Pdel[i].size();
659 MPI_Send(&outgoing, 1, MPI_INT, i, 0, MPI_COMM_WORLD);
660 while (Pdel[i].size() > 0) {
666 MPI_Send(&p_new, 1, pair_type, i, 0, MPI_COMM_WORLD);
671 while (Pdel[0].size() > 0) {
674 list<Critical_Pair_XPlor *>::iterator Ri = R.begin();
688 MPI_Recv(&incoming, 1, MPI_INT, 0, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
689 for (
int j = 0; j < incoming; ++j) {
690 MPI_Recv(&p_in, 1, pair_type, 0, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
691 list<Critical_Pair_XPlor *>::iterator Ri = R.begin();
694 ((*Ri)->first_index() != p_in.
first or
695 (*Ri)->second_index() != p_in.
second)
718 if (comm_id == 0) min_todo = P.size();
719 unsigned my_todo = R.size();
721 MPI_Reduce(&my_todo, &all_todo, 1, MPI_UNSIGNED, MPI_MAX, 0, MPI_COMM_WORLD);
723 min_todo = min_todo < all_todo ? all_todo : min_todo;
724 MPI_Bcast(&min_todo, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD);
737 MPI_Type_free(&pair_type);
738 cout << comm_id <<
" reduced " << number_of_spolys << endl;
739 MPI_Barrier(MPI_COMM_WORLD);
741 MPI_Reduce(&number_of_spolys, &total_spolys, 1, MPI_INT, MPI_SUM, 0, MPI_COMM_WORLD);
742 double max_bcast_time;
743 MPI_Reduce(&r_bcast_time, &max_bcast_time, 1, MPI_DOUBLE, MPI_MAX, 0, MPI_COMM_WORLD);
745 cout << total_spolys <<
" s-polynomials computed and reduced\n";
746 cout << max_bcast_time <<
" seconds spent in communication\n";
748 cout << G.size() <<
" polynomials before interreduction\n";
750 list<Abstract_Polynomial *> G_final;
752 G_final.push_back(g);
754 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
const COEF_TYPE * coefficient_array()
all the coefficients
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
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
const COEF_TYPE * numerator_array()
all the numerators
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
DEG_TYPE degree() const
returns the polynomial’s degree
Controls the creation of s-polynomials.
const UCOEF_TYPE * denominator_array()
all the denominators
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