PolyBoRi
minimal_elements.h
Go to the documentation of this file.
1 // -*- c++ -*-
2 //*****************************************************************************
14 //*****************************************************************************
15 
16 #ifndef polybori_groebner_minimal_elements_h_
17 #define polybori_groebner_minimal_elements_h_
18 
19 // include definitions
20 #include "groebner_defs.h"
21 #include "contained_variables.h"
22 //#include "nf.h"
23 
24 #include <vector>
25 #include <algorithm>
26 
28 MonomialSet mod_mon_set(const MonomialSet& as, const MonomialSet &vs); // nf.cc
29 MonomialSet mod_deg2_set(const MonomialSet& as, const MonomialSet &vs); // groebner_alg.cc
30 MonomialSet mod_var_set(const MonomialSet& as, const MonomialSet& vs); // groebner_alg.cc
31 
32 inline MonomialSet
34  if (s.isZero()) return s;
35  if (Polynomial(s).isOne()) return s;
37  int i=*nav;
38 
39 
40  if (Polynomial(s).hasConstantPart()) return MonomialSet(Polynomial(true, s.ring()));
41  int l=s.length();
42  if (l<=1) {
43  return s;
44  }
45 
46  if(l==2){
48  Monomial a=*it;
49  Monomial b=*(++it);
50  if (a.reducibleBy(b)) return MonomialSet(b.diagram());
51  else return s;
52  }
53 
54  MonomialSet s0_raw=s.subset0(i);
56  MonomialSet s1=minimal_elements_internal(s.subset1(i).diff(s0_raw));
57  if (!(s0.isZero())){
58  s1=s1.diff(s0.unateProduct(Polynomial(s1).usedVariablesExp().divisors(s.ring())));
59 
60  }
61  return s0.unite(s1.change(i));
62 
63 }
64 
65 inline MonomialSet
67  if (s.isZero()) return s;
68  if (Polynomial(s).isOne()) return s;
69 
70 
71 
72 
73  if (Polynomial(s).hasConstantPart())
74  return MonomialSet(Polynomial(true, s.ring()));
75  MonomialSet result(s.ring());
76  std::vector<idx_type> cv=contained_variables(s);
77  if ((cv.size()>0) && (s.length()==cv.size())){
78  return s;
79  } else {
80 
82  MonomialSet cv_set(s.ring());
83  for(z=cv.size()-1;z>=0;z--){
84  Monomial mv=Variable(cv[z], s.ring());
85  cv_set=cv_set.unite(mv.diagram());
86  }
87  for(z=0;z<cv.size();z++){
88  s=s.subset0(cv[z]);
89  }
90  result=cv_set;
91  }
92 
93  if (s.isZero()) return result;
94  PBORI_ASSERT(!(Polynomial(s).hasConstantPart()));
95 
96 
97 
99  idx_type i;
100  #if 1
101 
102  //first index of ZDD
103  i=*nav;
104  #else
105 
106  //first index of least lex. term
107  while(!(nav.isConstant())){
108  i=*nav;
109  nav.incrementElse();
110  }
111  #endif
112 
113 
114 
115  /*MonomialSet s0=minimal_elements_internal2(s.subset0(i));
116  MonomialSet s1=s.subset1(i);
117  if ((s0!=s1)&&(!(s1.diff(s0).isZero()))){
118  s1=minimal_elements_internal2(s1.unite(s0)).diff(s0);
119  } else return s0;
120  return s0.unite(s1.change(i));*/
121 
122  MonomialSet s0_raw=s.subset0(i);
124  MonomialSet s1=minimal_elements_internal2(s.subset1(i).diff(s0_raw));
125  if (!(s0.isZero())){
126  s1=s1.diff(s0.unateProduct(Polynomial(s1).usedVariables().divisors()));
127 
128  }
129  return s0.unite(s1.change(i)).unite(result);
130 
131 }
132 
133 inline std::vector<Exponent>
135  std::vector<Exponent> result;
136  if (s.isZero()) return result;
137  if ((Polynomial(s).isOne()) || (Polynomial(s).hasConstantPart())){
138  result.push_back(Exponent());
139  return result;
140  }
141  std::vector<idx_type> cv=contained_variables(s);
142  for(MonomialSet::size_type i=0;i<cv.size();i++){
143  s=s.subset0(cv[i]);
144  Exponent t;
145  t.insert(cv[i]);
146  result.push_back(t);
147  }
148 
149  if (s.isZero()){
150  return result;
151  } else {
152  std::vector<Exponent> exponents;
153  //Pol sp=s;
154  exponents.insert(exponents.end(), s.expBegin(),s.expEnd());
155  int nvars=s.ring().nVariables();
156  std::vector<std::vector<int> > occ_vecs(nvars);
157  for(MonomialSet::size_type i=0;i<exponents.size()-1;i++){
158  Exponent::const_iterator it=((const Exponent&) exponents[i]).begin();
159  Exponent::const_iterator end=((const Exponent&) exponents[i]).end();
160  while(it!=end){
161  occ_vecs[*it].push_back(i);
162  it++;
163  }
164  }
165  //now exponents are ordered
166  /*std::vector<std::set<int> > occ_sets(nvars);
167  for(i=occ_sets.size()-1;i>=0;i--){
168  occ_sets[i].insert(occ_vecs[i].begin(),occ_vecs[i].end());
169  }*/
170  std::vector<bool> still_minimal(exponents.size());
171  for(MonomialSet::size_type i=exponents.size()-1;i>=0;i--){
172  still_minimal[i]=true;
173  }
174 
175  //lex smalles is last so backwards
176  for(MonomialSet::size_type i=exponents.size()-1;i>=0;i--){
177  if (still_minimal[i]){
178  //we assume, that each exponents has deg>0
179  Exponent::const_iterator it=((const Exponent&) exponents[i]).begin();
180  Exponent::const_iterator end=((const Exponent&) exponents[i]).end();
181  int first_index=*it;
182  std::vector<int> occ_set=occ_vecs[first_index];
183  it++;
184  while(it!=end){
185 
186  std::vector<int> occ_set_next;
187  set_intersection(
188  occ_set.begin(),
189  occ_set.end(),
190  occ_vecs[*it].begin(),
191  occ_vecs[*it].end(),
192  std::back_insert_iterator<std::vector<int> >(occ_set_next));
193  occ_set=occ_set_next;
194  it++;
195  }
196 
197  std::vector<int>::const_iterator oc_it= occ_set.begin();
198  std::vector<int>::const_iterator oc_end= occ_set.end();
199  while(oc_it!=oc_end){
200  still_minimal[*oc_it]=false;
201  oc_it++;
202  }
203 
204  it=((const Exponent&) exponents[i]).begin();
205  while(it!=end){
206  std::vector<int> occ_set_difference;
207  set_difference(
208  occ_vecs[*it].begin(),
209  occ_vecs[*it].end(),
210  occ_set.begin(),
211  occ_set.end(),
212  std::back_insert_iterator<std::vector<int> >(occ_set_difference));
213  occ_vecs[*it]=occ_set_difference;
214  it++;
215  }
216  result.push_back(exponents[i]);
217  }
218  }
219 
220  }
221  return result;
222 }
223 
224 inline MonomialSet
226 #if 0
227  //return minimal_elements_internal2(s);
229 #else
230 #if 1
231  return s.minimalElements();
232 #else
234 
235  if (minElts!=s.minimalElements()){
236 
237  std::cout <<"ERROR minimalElements"<<std::endl;
238  std::cout <<"orig "<<s<< std::endl;
239  std::cout <<"correct "<<minElts<< std::endl;
240  std::cout <<"wrong "<<s.minimalElements()<< std::endl;
241  }
242 
243 
244  return minElts;
245 #endif
246 #endif
247 }
248 
249 
250 
251 inline MonomialSet
253 
254  if (m.isZero()) return m;
255 
256  if (m.ownsOne()) return Polynomial(1, m.ring()).diagram();
257 
261 
262 
263  typedef PBORI::CacheManager<CCacheTypes::minimal_elements>
264  cache_mgr_type;
265 
266 
267 
268  cache_mgr_type cache_mgr(m.ring());
269  PBORI::BoolePolynomial::navigator cached =
270  cache_mgr.find(m_nav);
271 
272 
273 
274  if (cached.isValid() ){
275  return cache_mgr.generate(cached);
276  }
277 
278  MonomialSet minimal_else=minimal_elements_cudd_style_unary(cache_mgr.generate(ms0));
279  MonomialSet minimal_then=minimal_elements_cudd_style_unary(mod_mon_set(cache_mgr.generate(ms1),minimal_else));
280  MonomialSet result(m.ring());
281  if ((minimal_else.navigation()==ms0) &&(minimal_then.navigation()==ms1)) result=m;
282  else
283  result= MonomialSet(*m_nav,minimal_then,minimal_else);//result0.unite(result1.change(index));
284 
285  cache_mgr.insert(m_nav, result.navigation());
286  return result;
287 }
288 
289 inline MonomialSet
291  Polynomial p_mod=mod;
292  if (m.isZero()) return m;
293  if (mod.ownsOne())
294  return MonomialSet(mod.ring());
295  if (m.ownsOne()) return Polynomial(1,m.ring()).diagram();
297  m=mod_var_set(m,mod_cv);
298  mod=mod_var_set(mod,mod_cv);
300  mod=mod_deg2_set(mod,mod_d2);
301  m=mod_deg2_set(m,mod_d2);
303  MonomialSet cv_orig=cv;
304  cv=cv.diff(mod);
305  mod=mod_var_set(mod,cv_orig);
306  m=mod_var_set(m,cv_orig);
307  m=m.diff(mod);
308  if (m.isZero()) return cv;
309  bool cv_empty=cv.isZero();
310 
311  MonomialSet result(m.ring());
312  int index=*m.navigation();
313 
314 
315 
316 
317  if (!mod.isZero())
318  {
319  MonomialSet::navigator nav_mod=mod.navigation();
320  while((!(nav_mod.isConstant())) && (index>*nav_mod)){
321  nav_mod.incrementElse();
322 
323  }
324  mod=MonomialSet(nav_mod, m.ring());
325  }
326 
330  MonomialSet::navigator mod_nav=mod.navigation();
331 
332  typedef PBORI::CacheManager<CCacheTypes::minimal_mod>
333  cache_mgr_type;
334 
335 
336 
337  cache_mgr_type cache_mgr(m.ring());
338  PBORI::BoolePolynomial::navigator cached =
339  cache_mgr.find(m_nav, mod_nav);
340 
341 
342 
343  if (cached.isValid() ){
344  return cv.unite((MonomialSet)cache_mgr.generate(cached));
345  }
346 
347  if (mod.isZero()){
348 
349  MonomialSet result0=do_minimal_elements_cudd_style(cache_mgr.generate(ms0),
350  mod);
351  MonomialSet result1= do_minimal_elements_cudd_style(cache_mgr.generate(ms1),
352  result0);
353  result= MonomialSet(index,result1,result0);//result0.unite(result1.change(index));
354 
355  } else {
356  if (*mod_nav>index){
358  result0=do_minimal_elements_cudd_style(cache_mgr.generate(ms0), mod);
360  cache_mgr.generate(ms1),result0.unite(mod));
361  if (result1.isZero()) {result=result0;}
362  else
363  {result= MonomialSet(index,result1,result0);}
364  } else {
365  PBORI_ASSERT(index==*mod_nav);
366  MonomialSet::navigator mod0=mod_nav.elseBranch();
367  MonomialSet::navigator mod1=mod_nav.thenBranch();
369  result0=do_minimal_elements_cudd_style(cache_mgr.generate(ms0), cache_mgr.generate(mod0));
370  //MonomialSet mod1=mod.subset1(index);
371  MonomialSet result1=
372  do_minimal_elements_cudd_style(cache_mgr.generate(ms1),
373  result0.unite(cache_mgr.generate(ms0).unite(cache_mgr.generate(mod1))));
374  result= MonomialSet(index,result1,result0);//result0.unite(result1.change(index));
375  }
376  }
377  cache_mgr.insert(m.navigation(), mod.navigation(), result.navigation());
378  if (cv_empty) return result;
379  else
380  return cv.unite(result);
381 }
382 
383 inline MonomialSet
386 }
387 
388 
389 
390 
391 
392 inline MonomialSet
394 
395  if (m.divisorsOf(lm).isZero()){
396  m = m.existAbstract(lm);
398 
399  return m.unateProduct(lm.set());
400  }
401  return lm.set() ;
402 }
403 
404 inline std::vector<Monomial>
406 
408  return std::vector<Monomial>(result.begin(), result.end());
409 }
410 
411 
412 
413 //#define MIN_ELEMENTS_BINARY 1
414 #ifdef MIN_ELEMENTS_BINARY
415 inline MonomialSet
417 
418  if (m.divisorsOf(lm).isZero()){
419  m=divide_monomial_divisors_out(m,lm);
420  //mod=divide_monomial_divisors_out(mod,lm);
421  return do_minimal_elements_cudd_style(m,mod);
422  }
423  return m.ring().one();
424 }
425 
426 
427 #else
428 
429 inline MonomialSet
431 
432  if (m.divisorsOf(lm).isZero()){
433 
434  m=m.existAbstract(lm);
435  mod=mod.existAbstract(lm);
436  //mod=divide_monomial_divisors_out(mod,lm);
437  m=mod_mon_set(m,mod);
438 
440  }
441  return m.ring().one();
442 }
443 #endif
444 
445 inline void
447  std::vector<Exponent>& result){
448  MonomialSet elements = minimal_elements_divided(m, lm, mod);
449  result.resize(elements.length());
450  std::copy(elements.expBegin(), elements.expEnd(), result.begin());
451 }
452 
453 
455 
456 #endif /* polybori_groebner_minimal_elements_h_ */
self minimalElements() const
Get minimal elements wrt. inclusion.
Definition: BooleSet.cc:207
MonomialSet contained_deg2_cudd_style(const MonomialSet &m)
Definition: contained_variables.h:63
dd_type one() const
Get decision diagram with all variables negated.
Definition: BoolePolyRing.cc:46
This class is just a wrapper for using variables for storing indices as interim data structure for Bo...
Definition: BooleExponent.h:34
MonomialSet do_minimal_elements_cudd_style(MonomialSet m, MonomialSet mod)
Definition: minimal_elements.h:290
std::vector< idx_type > contained_variables(const MonomialSet &m)
Definition: contained_variables.h:91
bool_type reducibleBy(const self &rhs) const
Test for reducibility.
Definition: BooleMonomial.h:190
bool_type isOne() const
Check whether polynomial is constant one.
Definition: BoolePolynomial.h:297
MonomialSet contained_variables_cudd_style(const MonomialSet &m)
Definition: contained_variables.h:25
polybori::BooleExponent Exponent
Definition: groebner_defs.h:33
set_type set() const
Get corresponding subset of of the powerset over all variables.
Definition: BooleMonomial.h:216
#define END_NAMESPACE_PBORIGB
Definition: groebner_defs.h:16
self existAbstract(const term_type &rhs) const
Definition: BooleSet.cc:252
MonomialSet minimal_elements_internal(const MonomialSet &s)
Definition: minimal_elements.h:33
bool_type hasConstantPart() const
Check whether polynomial has term one.
Definition: BoolePolynomial.h:303
self divisorsOf(const term_type &rhs) const
Compute intersection with divisors of rhs.
Definition: BooleSet.cc:156
MonomialSet mod_deg2_set(const MonomialSet &as, const MonomialSet &vs)
Definition: groebner_alg.cc:135
exp_iterator expEnd() const
Finish of iteration over exponent vectors.
Definition: BooleSet.cc:109
MonomialSet mod_var_set(const MonomialSet &as, const MonomialSet &vs)
Definition: groebner_alg.cc:85
BoolePolynomial Polynomial
Definition: embed.h:51
navigator navigation() const
Navigate through ZDD by incrementThen(), incrementElse(), and terminated()
Definition: CCuddDDFacade.h:455
#define BEGIN_NAMESPACE_PBORIGB
Definition: groebner_defs.h:15
self & insert(idx_type)
Insert variable with index idx in exponent vector.
Definition: BooleExponent.cc:171
self elseBranch() const
Increment in else direction.
Definition: CCuddNavigator.h:98
exp_iterator expBegin() const
Start of iteration over exponent vectors.
Definition: BooleSet.cc:101
This class wraps the underlying decicion diagram type and defines the necessary operations.
Definition: BoolePolynomial.h:85
const ring_type & ring() const
Get reference to ring.
Definition: CCuddDDFacade.h:250
polybori::BooleSet MonomialSet
Definition: groebner_defs.h:45
BooleVariable Variable
Definition: embed.h:52
bool isZero() const
Test whether diagram represents the empty set.
Definition: CCuddDDFacade.h:244
std::vector< Monomial > minimal_elements_multiplied(MonomialSet m, Monomial lm)
Definition: minimal_elements.h:405
MonomialSet minimal_elements_internal2(MonomialSet s)
Definition: minimal_elements.h:66
MonomialSet mod_mon_set(const MonomialSet &as, const MonomialSet &vs)
Definition: nf.cc:855
MonomialSet minimal_elements_cudd_style_unary(MonomialSet m)
Definition: minimal_elements.h:252
self thenBranch() const
Increment in then direction.
Definition: CCuddNavigator.h:92
#define PBORI_ASSERT(arg)
Definition: pbori_defs.h:118
MonomialSet minimal_elements(const MonomialSet &s)
Definition: minimal_elements.h:225
size_type nVariables() const
Get number of ring variables.
Definition: BoolePolyRing.h:126
bool_type ownsOne() const
Test whether the empty set is included.
Definition: BooleSet.h:211
MonomialSet minimal_elements_cudd_style(MonomialSet m)
Definition: minimal_elements.h:384
const dd_type & diagram() const
Read-only access to internal decision diagramm structure.
Definition: BooleMonomial.h:213
const_iterator end() const
Finish of iteration over terms.
Definition: BooleSet.cc:78
Definition: BoolePolynomial.h:70
self change(idx_type idx) const
Definition: BooleSet.h:169
void minimal_elements_divided(MonomialSet m, Monomial lm, MonomialSet mod, std::vector< Exponent > &result)
Definition: minimal_elements.h:446
std::vector< Exponent > minimal_elements_internal3(MonomialSet s)
Definition: minimal_elements.h:134
polybori::CTypes::idx_type idx_type
Definition: groebner_defs.h:44
std::size_t size_type
Type for lengths, dimensions, etc.
Definition: pbori_defs.h:219
MonomialSet minimal_elements_multiplied_monoms(MonomialSet m, Monomial lm)
Definition: minimal_elements.h:393
data_type::const_iterator const_iterator
Definition: BooleExponent.h:52
bool_type isConstant() const
Check whether constant node was reached.
Definition: CCuddNavigator.h:172
This class defines an iterator for navigating through then and else branches of ZDDs.
Definition: CCuddNavigator.h:36
self & incrementElse()
Increment in else direction.
Definition: CCuddNavigator.h:203
const_iterator begin() const
Start of iteration over terms.
Definition: BooleSet.cc:70
size_type length() const
Returns number of terms (deprecated)
Definition: BooleSet.h:245
Definition: BooleSet.h:57
This class is just a wrapper for using variables from cudd's decicion diagram.
Definition: BooleMonomial.h:50
self & push_back(idx_type idx)
Insert variable with index idx in exponent vector (trying end first)
Definition: BooleExponent.cc:188