PolyBoRi
pbori_routines_misc.h
Go to the documentation of this file.
00001 // -*- c++ -*-
00002 //*****************************************************************************
00013 //*****************************************************************************
00014 
00015 #ifndef polybori_routines_pbori_routines_misc_h_
00016 #define polybori_routines_pbori_routines_misc_h_
00017 
00018 // include basic definitions
00019 #include <polybori/pbori_defs.h>
00020 
00021 // temprarily
00022 #include <polybori/cache/CacheManager.h>
00023 
00024 #include <polybori/diagram/CDDOperations.h>
00025 
00026 BEGIN_NAMESPACE_PBORI
00027 
00028 template<class Iterator>
00029 typename Iterator::value_type
00030 index_vector_hash(Iterator start, Iterator finish){
00031 
00032   typedef typename Iterator::value_type value_type;
00033 
00034   value_type vars = 0;
00035   value_type sum = 0;
00036  
00037   while (start != finish){
00038     vars++;
00039     sum += ((*start)+1) * ((*start)+1);
00040     ++start;
00041   }
00042   return sum * vars;
00043 }
00044 
00047 template <class DegreeCacher, class NaviType>
00048 typename NaviType::deg_type
00049 dd_cached_degree(const DegreeCacher& cache, NaviType navi) {
00050 
00051   typedef typename NaviType::deg_type deg_type;
00052 
00053   if (navi.isConstant()){ // No need for caching of constant nodes' degrees
00054     if (navi.terminalValue())
00055         return 0;
00056     else
00057         return -1;
00058   }
00059  
00060   // Look whether result was cached before
00061   typename DegreeCacher::node_type result = cache.find(navi);
00062   if (result.isValid())
00063     return (*result);
00064 
00065   // Get degree of then branch (contains at least one valid path)...
00066   deg_type deg = dd_cached_degree(cache, navi.thenBranch()) + 1;
00067  
00068   // ... combine with degree of else branch
00069   deg = std::max(deg,  dd_cached_degree(cache, navi.elseBranch()) );
00070 
00071   // Write result to cache
00072   cache.insert(navi, deg);
00073  
00074   return deg;
00075 }
00076 
00081 template <class DegreeCacher, class NaviType, class SizeType>
00082 typename NaviType::deg_type
00083 dd_cached_degree(const DegreeCacher& cache, NaviType navi, SizeType bound) {
00084 
00085   typedef typename NaviType::deg_type deg_type;
00086 
00087   if (bound < 0) {
00088     return -1;
00089   }
00090 
00091   // No need for caching of constant nodes' degrees
00092   if (bound == 0) {
00093     while(!navi.isConstant())
00094       navi.incrementElse();
00095   }
00096   if (navi.isConstant())
00097     return (navi.terminalValue()? 0: -1);
00098 
00099   // Look whether result was cached before
00100   typename DegreeCacher::node_type result = cache.find(navi, bound);
00101   if (result.isValid())
00102     return *result;
00103 
00104   // Get degree of then branch (contains at least one valid path)...
00105   deg_type deg = dd_cached_degree(cache, navi.thenBranch(), bound - 1);
00106   if (deg >= 0) ++deg;
00107 
00108   // ... combine with degree of else branch
00109   if (bound > deg)              // if deg <= bound, we are already finished
00110     deg = std::max(deg,  dd_cached_degree(cache, navi.elseBranch(), bound) );
00111 
00112   // Write result to cache
00113   if (deg >= 0) {
00114     PBORI_ASSERT(bound >= 0);
00115     cache.insert(navi, bound, deg);
00116   }
00117 
00118   return deg;
00119 }
00120 
00121 template <class Iterator, class NameGenerator, 
00122           class Separator, class EmptySetType, 
00123           class OStreamType>
00124 void 
00125 dd_print_term(Iterator start, Iterator finish, const NameGenerator& get_name,
00126               const Separator& sep, const EmptySetType& emptyset, 
00127               OStreamType& os){
00128 
00129   if (start != finish){
00130     os << get_name(*start);
00131     ++start;
00132   }
00133   else
00134     os << emptyset();
00135 
00136   while (start != finish){
00137     os << sep() << get_name(*start);
00138     ++start;
00139   }
00140 }
00141 
00142 template <class TermType, class NameGenerator,
00143           class Separator, class EmptySetType,
00144           class OStreamType>
00145 void 
00146 dd_print_term(const TermType& term, const NameGenerator& get_name,
00147               const Separator& sep, const EmptySetType& emptyset, 
00148               OStreamType& os){
00149   dd_print_term(term.begin(), term.end(), get_name, sep, emptyset, os);
00150 }
00151 
00152 
00153 template <class Iterator, class NameGenerator,
00154           class Separator, class InnerSeparator, 
00155           class EmptySetType, class OStreamType>
00156 void 
00157 dd_print_terms(Iterator start, Iterator finish, const NameGenerator& get_name,
00158                const Separator& sep, const InnerSeparator& innersep,
00159                const EmptySetType& emptyset, OStreamType& os) {
00160 
00161   if (start != finish){
00162     dd_print_term(*start, get_name, innersep, emptyset, os);
00163     ++start;
00164   }
00165 
00166   while (start != finish){
00167     os << sep(); 
00168     dd_print_term(*start, get_name, innersep, emptyset, os);
00169     ++start;
00170   }
00171 
00172 }
00173 
00174 
00175 template <bool use_fast, 
00176           class CacheType, class NaviType, class PolyType>
00177 PolyType
00178 dd_multiply(const CacheType& cache_mgr,
00179             NaviType firstNavi, NaviType secondNavi, PolyType init){
00180 
00181   // Extract subtypes
00182   typedef typename PolyType::dd_type dd_type;
00183   typedef typename NaviType::idx_type idx_type;
00184   typedef NaviType navigator;
00185 
00186   if (firstNavi.isConstant()) {
00187     if(firstNavi.terminalValue())
00188       return cache_mgr.generate(secondNavi);
00189     else 
00190       return cache_mgr.zero();
00191   }
00192 
00193   if (secondNavi.isConstant()) {
00194     if(secondNavi.terminalValue())
00195       return cache_mgr.generate(firstNavi);
00196     else 
00197       return cache_mgr.zero();
00198   }
00199   if (firstNavi == secondNavi)
00200     return cache_mgr.generate(firstNavi);
00201   
00202   // Look up, whether operation was already used
00203   navigator cached = cache_mgr.find(firstNavi, secondNavi);
00204   PolyType result(cache_mgr.zero());
00205 
00206   if (cached.isValid()) {       // Cache lookup sucessful
00207     return cache_mgr.generate(cached);
00208   }
00209   else {                        // Cache lookup not sucessful
00210     // Get top variable's index
00211 
00212     if (*secondNavi < *firstNavi)
00213       std::swap(firstNavi, secondNavi);
00214 
00215     idx_type index = *firstNavi;
00216 
00217     // Get then- and else-branches wrt. current indexed variable
00218     navigator as0 = firstNavi.elseBranch();
00219     navigator as1 = firstNavi.thenBranch();
00220 
00221     navigator bs0;
00222     navigator bs1;
00223 
00224     if (*secondNavi == index) {
00225       bs0 = secondNavi.elseBranch();
00226       bs1 = secondNavi.thenBranch();
00227     }
00228     else {
00229       bs0 = secondNavi;
00230       bs1 = cache_mgr.zero().navigation();
00231     }
00232     PolyType result0 = dd_multiply<use_fast>(cache_mgr, as0, bs0, init);
00233     PolyType result1(cache_mgr.zero());
00234 
00235     // use fast multiplication
00236     if (use_fast && (*firstNavi == *secondNavi)) {
00237 
00238       PolyType res10 = PolyType(cache_mgr.generate(as1)) +
00239         PolyType(cache_mgr.generate(as0));
00240       PolyType res01 = PolyType(cache_mgr.generate(bs0)) +
00241         PolyType(cache_mgr.generate(bs1));
00242       
00243       result1 = dd_multiply<use_fast>(cache_mgr, res10.navigation(),
00244                                         res01.navigation(), init) - result0;
00245     } 
00246     // not using fast multiplication
00247     else if (as0 == as1) {
00248       result1 = dd_multiply<use_fast>(cache_mgr, bs0, as1, init);
00249 
00250     }
00251     else {
00252       result1 = dd_multiply<use_fast>(cache_mgr, as0, bs1, init); 
00253       if (bs0 != bs1){
00254         PolyType bs01 = PolyType(cache_mgr.generate(bs0)) + 
00255           PolyType(cache_mgr.generate(bs1));
00256         
00257         result1 += 
00258           dd_multiply<use_fast>(cache_mgr, bs01.navigation(), as1, init);
00259       }
00260     }
00261     result = dd_type(index, result1.diagram(), result0.diagram());
00262 
00263     // Insert in cache
00264     cache_mgr.insert(firstNavi, secondNavi, result.navigation());
00265   } // end of Cache lookup not sucessful
00266 
00267   return result;
00268 }
00269 
00270 template <class CacheType, class NaviType, class PolyType>
00271 PolyType
00272 dd_multiply_recursively(const CacheType& cache_mgr,
00273                         NaviType firstNavi, NaviType secondNavi, PolyType init){
00274 
00275   enum { use_fast = 
00276 #ifdef PBORI_FAST_MULTIPLICATION
00277          true
00278 #else
00279          false
00280 #endif
00281   };
00282 
00283   return dd_multiply<use_fast>(cache_mgr, firstNavi, secondNavi, init);
00284 }
00285 
00286 template <class CacheType, class NaviType, class PolyType>
00287 PolyType
00288 dd_multiply_recursively_monom(const CacheType& cache_mgr,
00289                         NaviType monomNavi, NaviType navi, PolyType init){
00290 
00291   // Extract subtypes
00292   typedef typename PolyType::dd_type dd_type;
00293   typedef typename NaviType::idx_type idx_type;
00294   typedef NaviType navigator;
00295 
00296   if (monomNavi.isConstant()) {
00297     if(monomNavi.terminalValue())
00298       return cache_mgr.generate(navi);
00299     else 
00300       return cache_mgr.zero();
00301   }
00302 
00303   PBORI_ASSERT(monomNavi.elseBranch().isEmpty());
00304 
00305   if (navi.isConstant()) {
00306     if(navi.terminalValue())
00307       return cache_mgr.generate(monomNavi);
00308     else 
00309       return cache_mgr.zero();
00310   }
00311   if (monomNavi == navi)
00312     return cache_mgr.generate(monomNavi);
00313   
00314   // Look up, whether operation was already used
00315   navigator cached = cache_mgr.find(monomNavi, navi);
00316 
00317   if (cached.isValid()) {       // Cache lookup sucessful
00318     return cache_mgr.generate(cached);
00319   }
00320 
00321   // Cache lookup not sucessful
00322   // Get top variables' index
00323 
00324   idx_type index = *navi;
00325   idx_type monomIndex = *monomNavi;
00326 
00327   if (monomIndex < index) {     // Case: index may occure within monom
00328     init = dd_multiply_recursively_monom(cache_mgr, monomNavi.thenBranch(), navi,
00329                                    init).diagram().change(monomIndex);
00330   }
00331   else if (monomIndex == index) { // Case: monom and poly start with same index
00332 
00333     // Increment navigators
00334     navigator monomThen = monomNavi.thenBranch();
00335     navigator naviThen = navi.thenBranch();
00336     navigator naviElse = navi.elseBranch();
00337 
00338     if (naviThen != naviElse)
00339       init = (dd_multiply_recursively_monom(cache_mgr, monomThen, naviThen, init)
00340               + dd_multiply_recursively_monom(cache_mgr, monomThen, naviElse,
00341                                               init)).diagram().change(index);
00342   }
00343   else {                        // Case: var(index) not part of monomial
00344     
00345     init = 
00346       dd_type(index,  
00347               dd_multiply_recursively_monom(cache_mgr, monomNavi, navi.thenBranch(), 
00348                                       init).diagram(),
00349               dd_multiply_recursively_monom(cache_mgr, monomNavi, navi.elseBranch(),
00350                                       init).diagram() );
00351   }
00352   
00353   // Insert in cache
00354   cache_mgr.insert(monomNavi, navi, init.navigation());
00355 
00356   return init;
00357 }
00358 
00359 
00360 template <class DDGenerator, class Iterator, class NaviType, class PolyType>
00361 PolyType
00362 dd_multiply_recursively_exp(const DDGenerator& ddgen,
00363                             Iterator start, Iterator finish,
00364                             NaviType navi, PolyType init){
00365   // Extract subtypes
00366   typedef typename NaviType::idx_type idx_type;
00367   typedef typename PolyType::dd_type dd_type;
00368   typedef NaviType navigator;
00369 
00370   if (start == finish)
00371     return ddgen.generate(navi);
00372 
00373   PolyType result(ddgen.zero());
00374   if (navi.isConstant()) {
00375     if(navi.terminalValue()) {
00376 
00377       std::reverse_iterator<Iterator> rstart(finish), rfinish(start);
00378       result = ddgen.generate(navi);
00379       while (rstart != rfinish) {
00380         result = result.diagram().change(*rstart);
00381         ++rstart;
00382       }
00383     }
00384     else 
00385       return ddgen.zero();
00386 
00387     return result;
00388   }
00389 
00390   // Cache lookup not sucessful
00391   // Get top variables' index
00392 
00393   idx_type index = *navi;
00394   idx_type monomIndex = *start;
00395 
00396   if (monomIndex < index) {     // Case: index may occure within monom
00397 
00398     Iterator next(start);
00399     while( (next != finish) && (*next < index) )
00400       ++next;
00401 
00402     result = dd_multiply_recursively_exp(ddgen, next, finish, navi, init);
00403 
00404     std::reverse_iterator<Iterator> rstart(next), rfinish(start);
00405     while (rstart != rfinish) {
00406       result = result.diagram().change(*rstart);
00407       ++rstart;
00408     }
00409   }
00410   else if (monomIndex == index) { // Case: monom and poly start with same index
00411 
00412     // Increment navigators
00413     ++start;
00414 
00415     navigator naviThen = navi.thenBranch();
00416     navigator naviElse = navi.elseBranch();
00417 
00418     if (naviThen != naviElse)
00419       result =(dd_multiply_recursively_exp(ddgen, start, finish, naviThen, init)
00420               + dd_multiply_recursively_exp(ddgen, start, finish, naviElse,
00421                                             init)).diagram().change(index);
00422   }
00423   else {                        // Case: var(index) not part of monomial
00424     
00425     result = 
00426       dd_type(index,  
00427               dd_multiply_recursively_exp(ddgen, start, finish,
00428                                           navi.thenBranch(), init).diagram(),
00429               dd_multiply_recursively_exp(ddgen, start, finish, 
00430                                           navi.elseBranch(), init).diagram() );
00431   }
00432 
00433   return result;
00434 }
00435 
00436 template<class DegCacheMgr, class NaviType, class SizeType>
00437 bool max_degree_on_then(const DegCacheMgr& deg_mgr, NaviType navi,
00438                         SizeType degree, valid_tag is_descending) {
00439 
00440   navi.incrementThen();
00441   return ((dd_cached_degree(deg_mgr, navi, degree - 1) + 1) == degree);
00442 }
00443 
00444 template<class DegCacheMgr, class NaviType, class SizeType>
00445 bool max_degree_on_then(const DegCacheMgr& deg_mgr, NaviType navi,
00446                         SizeType degree, invalid_tag non_descending) {
00447 
00448   navi.incrementElse();
00449   return (dd_cached_degree(deg_mgr, navi, degree) != degree);
00450 }
00451 
00452 template <class NaviType>
00453 NaviType dd_get_constant(NaviType navi) {
00454   while(!navi.isConstant()) {
00455     navi.incrementElse();
00456   }
00457 
00458   return navi;
00459 }
00460 
00461 template <class T> class CIndexCacheHandle;
00462 // with degree bound
00463 template <class CacheType, class DegCacheMgr, class NaviType, 
00464           class TermType, class SizeType, class DescendingProperty>
00465 TermType
00466 dd_recursive_degree_lead(const CacheType& cache_mgr, const DegCacheMgr&
00467                          deg_mgr,
00468                          NaviType navi, TermType init, SizeType degree,
00469                          DescendingProperty prop) {
00470 
00471 
00472   typedef CIndexCacheHandle<NaviType> deg2node;
00473 
00474   if PBORI_UNLIKELY(degree < 0)
00475     return cache_mgr.zero();
00476 
00477   if (navi.isConstant())
00478     return cache_mgr.generate(navi);
00479 
00480   if (degree == 0)
00481     return cache_mgr.generate(dd_get_constant(navi));
00482 
00483   // Check cache for previous results
00484   NaviType cached = cache_mgr.find(navi, deg2node(degree, cache_mgr.manager()));
00485   if (cached.isValid())
00486     return cache_mgr.generate(cached);
00487 
00488   // Go to next branch
00489   if ( max_degree_on_then(deg_mgr, navi, degree, prop) ) {
00490     NaviType then_branch = navi.thenBranch();
00491     init = dd_recursive_degree_lead(cache_mgr, deg_mgr, then_branch, 
00492         init, degree - 1, prop);
00493 
00494     if  ((navi.elseBranch().isEmpty()) && (init.navigation()==then_branch))
00495       init = cache_mgr.generate(navi);
00496     else
00497       init = init.change(*navi);
00498                                     
00499   }
00500   else {
00501     init = dd_recursive_degree_lead(cache_mgr, deg_mgr, navi.elseBranch(), 
00502                                     init, degree, prop);
00503   }
00505   NaviType resultNavi(init.navigation());
00506   cache_mgr.insert(navi, deg2node(degree, cache_mgr.manager()), resultNavi);
00507 //   if (!cache_mgr.find(resultNavi).isValid())
00508 //     deg_mgr.insert(resultNavi, degree, degree);
00509 
00510   return init;
00511 }
00512 
00513 template <class CacheType, class DegCacheMgr, class NaviType, 
00514           class TermType, class DescendingProperty>
00515 TermType
00516 dd_recursive_degree_lead(const CacheType& cache_mgr, const DegCacheMgr& deg_mgr,
00517                          NaviType navi, TermType init, DescendingProperty prop){
00518 
00519 //   if (navi.isConstant())
00520 //     return cache_mgr.generate(dd_get_constant(navi));
00521   
00522   return dd_recursive_degree_lead(cache_mgr, deg_mgr, navi, init,
00523                                   dd_cached_degree(deg_mgr, navi), prop);
00524 }
00525 
00526 template <class CacheType, class DegCacheMgr, class NaviType, 
00527           class TermType, class SizeType, class DescendingProperty>
00528 TermType&
00529 dd_recursive_degree_leadexp(const CacheType& cache_mgr, 
00530                             const DegCacheMgr& deg_mgr,
00531                             NaviType navi, TermType& result,  
00532                             SizeType degree,
00533                             DescendingProperty prop) {
00534   typedef CIndexCacheHandle<NaviType> deg2node;
00535   if PBORI_UNLIKELY(degree < 0) {
00536     result.resize(0);
00537     return result;
00538   }
00539   if (navi.isConstant())
00540     return result;
00541 
00542   if (degree == 0) {
00543     while (!navi.isConstant())
00544       navi.incrementElse();
00545     if (!navi.terminalValue())
00546       result.resize(0);
00547 
00548     return result;
00549   }
00550 
00551   // Check cache for previous result
00552   NaviType cached = cache_mgr.find(navi, deg2node(degree, cache_mgr.manager()));
00553   if (cached.isValid())
00554     return result = result.multiplyFirst(cache_mgr.generate(cached));
00555 
00556   // Prepare next branch
00557   if ( max_degree_on_then(deg_mgr, navi, degree, prop) ) {
00558     result.push_back(*navi);
00559     navi.incrementThen();
00560     --degree;
00561   }
00562   else 
00563     navi.incrementElse();
00564 
00565   return 
00566     dd_recursive_degree_leadexp(cache_mgr, deg_mgr, navi, result, degree, prop);
00567 }
00568 
00569 template <class CacheType, class DegCacheMgr, class NaviType, 
00570           class TermType, class DescendingProperty>
00571 TermType&
00572 dd_recursive_degree_leadexp(const CacheType& cache_mgr, 
00573                             const DegCacheMgr& deg_mgr,
00574                             NaviType navi, TermType& result,  
00575                             DescendingProperty prop) {
00576 
00577   if (navi.isConstant())
00578     return result;
00579 
00580   return dd_recursive_degree_leadexp(cache_mgr, deg_mgr, navi, result, 
00581                                      dd_cached_degree(deg_mgr, navi), prop);
00582 }
00583 
00584 // Existential Abstraction. Given a ZDD F, and a set of variables
00585 // S, remove all occurrences of s in S from any subset in F. This can
00586 // be implemented by cofactoring F with respect to s = 1 and s = 0,
00587 // then forming the union of these results. 
00588 
00589 
00590 template <class CacheType, class NaviType, class TermType>
00591 TermType
00592 dd_existential_abstraction(const CacheType& cache_mgr, 
00593                            NaviType varsNavi, NaviType navi, TermType init){
00594 
00595   typedef typename TermType::dd_type dd_type;
00596   typedef typename TermType::idx_type idx_type;
00597 
00598   if (navi.isConstant()) 
00599     return cache_mgr.generate(navi);
00600 
00601   idx_type index(*navi);
00602   while (!varsNavi.isConstant() && ((*varsNavi) < index))
00603     varsNavi.incrementThen();
00604 
00605   if (varsNavi.isConstant())
00606     return cache_mgr.generate(navi);
00607 
00608   // Check cache for previous result
00609   NaviType cached = cache_mgr.find(varsNavi, navi);
00610   if (cached.isValid()) 
00611     return cache_mgr.generate(cached);
00612 
00613   NaviType thenNavi(navi.thenBranch()), elseNavi(navi.elseBranch()); 
00614 
00615   TermType thenResult = 
00616     dd_existential_abstraction(cache_mgr, varsNavi, thenNavi, init);
00617 
00618   TermType elseResult = 
00619     dd_existential_abstraction(cache_mgr, varsNavi, elseNavi, init);
00620 
00621   if ((*varsNavi) == index)
00622     init = thenResult.unite(elseResult);
00623   else if ( (thenResult.navigation() == thenNavi) && 
00624             (elseResult.navigation() == elseNavi)  )
00625     init = cache_mgr.generate(navi);
00626   else
00627     init = dd_type(index, thenResult, elseResult);
00628 
00629   // Insert result to cache
00630   cache_mgr.insert(varsNavi, navi, init.navigation());
00631 
00632   return init;
00633 }
00634 
00635 
00636 
00637 template <class CacheType, class NaviType, class PolyType>
00638 PolyType
00639 dd_divide_recursively(const CacheType& cache_mgr,
00640                       NaviType navi, NaviType monomNavi, PolyType init){
00641 
00642   // Extract subtypes
00643   typedef typename NaviType::idx_type idx_type;
00644   typedef NaviType navigator;
00645   typedef typename PolyType::dd_type dd_type;
00646 
00647   if (monomNavi.isConstant()) {
00648     PBORI_ASSERT(monomNavi.terminalValue() == true);
00649     return cache_mgr.generate(navi);
00650   }
00651 
00652   PBORI_ASSERT(monomNavi.elseBranch().isEmpty());
00653 
00654   if (navi.isConstant()) 
00655     return cache_mgr.zero();
00656 
00657   if (monomNavi == navi)
00658     return cache_mgr.one();
00659   
00660   // Look up, whether operation was already used
00661   navigator cached = cache_mgr.find(navi, monomNavi);
00662 
00663   if (cached.isValid()) {       // Cache lookup sucessful
00664     return cache_mgr.generate(cached);
00665   }
00666 
00667   // Cache lookup not sucessful
00668   // Get top variables' index
00669 
00670   idx_type index = *navi;
00671   idx_type monomIndex = *monomNavi;
00672 
00673   if (monomIndex == index) {    // Case: monom and poly start with same index
00674 
00675     // Increment navigators
00676     navigator monomThen =  monomNavi.thenBranch();
00677     navigator naviThen = navi.thenBranch();
00678 
00679     init = dd_divide_recursively(cache_mgr, naviThen, monomThen, init);
00680   }
00681   else if (monomIndex > index) { // Case: monomIndex may occure within poly
00682     
00683     init = 
00684       dd_type(index,  
00685               dd_divide_recursively(cache_mgr,  navi.thenBranch(), monomNavi,
00686                                       init).diagram(),
00687               dd_divide_recursively(cache_mgr, navi.elseBranch(), monomNavi, 
00688                                       init).diagram() );
00689   }
00690   
00691   // Insert in cache
00692   cache_mgr.insert(navi, monomNavi,  init.navigation());
00693 
00694   return init;
00695 }
00696 
00697 
00698 
00699 template <class DDGenerator, class Iterator, class NaviType, class PolyType>
00700 PolyType
00701 dd_divide_recursively_exp(const DDGenerator& ddgen,
00702                           NaviType navi, Iterator start, Iterator finish,
00703                           PolyType init){
00704 
00705   // Extract subtypes
00706   typedef typename NaviType::idx_type idx_type;
00707   typedef typename PolyType::dd_type dd_type;
00708   typedef NaviType navigator;
00709 
00710   if (start == finish)
00711     return ddgen.generate(navi);
00712 
00713   if (navi.isConstant()) 
00714     return ddgen.zero();
00715   
00716 
00717   // Cache lookup not sucessful
00718   // Get top variables' index
00719 
00720   idx_type index = *navi;
00721   idx_type monomIndex = *start;
00722 
00723   PolyType result(ddgen.zero());
00724   if (monomIndex == index) {    // Case: monom and poly start with same index
00725 
00726     // Increment navigators
00727     ++start;
00728     navigator naviThen = navi.thenBranch();
00729 
00730     result = dd_divide_recursively_exp(ddgen, naviThen, start, finish, init);
00731   }
00732   else if (monomIndex > index) { // Case: monomIndex may occure within poly
00733     
00734     result = 
00735       dd_type(index,  
00736               dd_divide_recursively_exp(ddgen, navi.thenBranch(), start, finish,
00737                                         init).diagram(),
00738               dd_divide_recursively_exp(ddgen, navi.elseBranch(), start, finish,
00739                                         init).diagram() );
00740   }
00741   else 
00742     result = ddgen.zero();
00743   
00744   return result;
00745 }
00746 
00749 template <class CacheType, class NaviType, class MonomType>
00750 MonomType
00751 cached_used_vars(const CacheType& cache, NaviType navi, MonomType init) {
00752 
00753   if (navi.isConstant()) // No need for caching of constant nodes' degrees
00754     return init;
00755  
00756   // Look whether result was cached before
00757   NaviType cached_result = cache.find(navi);
00758 
00759   typedef typename MonomType::poly_type poly_type;
00760   if (cached_result.isValid())
00761     return CDDOperations<typename MonomType::dd_type, 
00762     MonomType>().getMonomial(cache.generate(cached_result));
00763   
00764   MonomType result = cached_used_vars(cache, navi.thenBranch(), init);
00765   result *= cached_used_vars(cache, navi.elseBranch(), init);
00766 
00767   result = result.change(*navi);
00768 
00769   // Write result to cache
00770   cache.insert(navi, result.diagram().navigation());
00771  
00772   return result;
00773 }
00774 
00775 template <class NaviType, class Iterator>
00776 bool
00777 dd_owns(NaviType navi, Iterator start, Iterator finish) {
00778 
00779   if (start == finish) {
00780     while(!navi.isConstant())
00781       navi.incrementElse();
00782     return navi.terminalValue();
00783   }
00784 
00785   while(!navi.isConstant() && (*start > *navi) )
00786     navi.incrementElse();
00787 
00788   if (navi.isConstant() || (*start != *navi))
00789     return false;
00790 
00791   return dd_owns(navi.thenBranch(), ++start, finish);
00792 }
00793 
00796 template <class NaviType, class MonomIterator>
00797 bool
00798 dd_contains_divs_of_dec_deg(NaviType navi, 
00799                             MonomIterator start, MonomIterator finish) {
00800 
00801   // Managing trivial cases
00802 
00803   if (start == finish)          // divisors of monomial 1 is the empty set,
00804                                 // which is always contained
00805     return true;
00806 
00807   if (navi.isConstant()) {      // set is empty or owns 1 only
00808     if (navi.terminalValue())   // set == {1}
00809       return (++start == finish);  // whether monom is of degree one
00810     else                        // empty set contains no divisors
00811       return false;
00812 
00813   }
00814 
00815   // Actual computations
00816 
00817   if (*navi < *start)
00818     return dd_contains_divs_of_dec_deg(navi.elseBranch(), start, finish);
00819 
00820 
00821   if (*navi > *start) {
00822     if (++start == finish)      // if monom is of degree one
00823       return owns_one(navi);
00824     else                        
00825       return false;
00826   }  
00827 
00828   // (*navi == *start) here
00829   ++start;
00830   return dd_owns(navi.elseBranch(), start, finish) && 
00831     dd_contains_divs_of_dec_deg(navi.thenBranch(), start, finish);
00832 }
00833 
00834 
00835 
00836 // determine the part of a polynomials of a given degree
00837 template <class CacheType, class NaviType, class DegType, class SetType>
00838 SetType
00839 dd_graded_part(const CacheType& cache, NaviType navi, DegType deg,  
00840                SetType init) {
00841 
00842 
00843   if (deg == 0) {
00844     while(!navi.isConstant())
00845       navi.incrementElse();
00846     return cache.generate(navi);
00847   }
00848 
00849   if(navi.isConstant())
00850     return cache.zero();
00851 
00852   // Look whether result was cached before
00853   NaviType cached = cache.find(navi, deg);
00854 
00855   if (cached.isValid())
00856     return cache.generate(cached);
00857 
00858   SetType result = 
00859     SetType(*navi,  
00860             dd_graded_part(cache, navi.thenBranch(), deg - 1, init),
00861             dd_graded_part(cache, navi.elseBranch(), deg, init)
00862             );
00863 
00864   // store result for later reuse
00865   cache.insert(navi, deg, result.navigation());
00866 
00867   return result;
00868 }
00869 
00870 
00874 template <class CacheManager, class NaviType, class SetType>
00875 SetType
00876 dd_first_divisors_of(CacheManager cache_mgr, NaviType navi, 
00877                      NaviType rhsNavi, SetType init ) {
00878 
00879   typedef typename SetType::dd_type dd_type;
00880   while( (!navi.isConstant()) && (*rhsNavi != *navi) ) {
00881 
00882     if ( (*rhsNavi < *navi) && (!rhsNavi.isConstant()) ) 
00883       rhsNavi.incrementThen();
00884     else 
00885       navi.incrementElse();  
00886   }
00887 
00888   if (navi.isConstant())        // At end of path
00889     return cache_mgr.generate(navi);
00890 
00891   // Look up, whether operation was already used
00892   NaviType result = cache_mgr.find(navi, rhsNavi);
00893     
00894   if (result.isValid())       // Cache lookup sucessful
00895     return  cache_mgr.generate(result);
00896   
00897   PBORI_ASSERT(*rhsNavi == *navi);
00898    
00899   // Compute new result
00900   init = dd_type(*rhsNavi,  
00901                  dd_first_divisors_of(cache_mgr, navi.thenBranch(), rhsNavi, 
00902                                       init).diagram(),
00903                  dd_first_divisors_of(cache_mgr, navi.elseBranch(), rhsNavi, 
00904                                       init).diagram() );
00905   // Insert result to cache
00906   cache_mgr.insert(navi, rhsNavi, init.navigation());
00907 
00908   return init;
00909 }
00910 
00911 template <class CacheType, class NaviType, class SetType>
00912 SetType
00913 dd_first_multiples_of(const CacheType& cache_mgr,
00914                       NaviType navi, NaviType rhsNavi, SetType init){
00915 
00916   typedef typename SetType::dd_type dd_type;
00917 
00918   if(rhsNavi.isConstant()) {
00919     PBORI_ASSERT(rhsNavi.terminalValue() == true);
00920     return cache_mgr.generate(navi);
00921   }
00922 
00923   if (navi.isConstant() || (*navi > *rhsNavi)) 
00924     return cache_mgr.zero();
00925 
00926   if (*navi == *rhsNavi)
00927     return dd_first_multiples_of(cache_mgr, navi.thenBranch(), 
00928                                  rhsNavi.thenBranch(), init).change(*navi);
00929 
00930   // Look up old result - if any
00931   NaviType result = cache_mgr.find(navi, rhsNavi);
00932 
00933   if (result.isValid())
00934     return cache_mgr.generate(result);
00935 
00936   // Compute new result
00937   init = dd_type(*navi,
00938                   dd_first_multiples_of(cache_mgr, navi.thenBranch(), 
00939                                         rhsNavi, init).diagram(),
00940                   dd_first_multiples_of(cache_mgr, navi.elseBranch(), 
00941                                         rhsNavi, init).diagram() );
00942 
00943   // Insert new result in cache
00944   cache_mgr.insert(navi, rhsNavi, init.navigation());
00945 
00946   return init;
00947 }
00948 
00949 
00952 template <class PolyType, class RingType, class MapType, class NaviType>
00953 PolyType
00954 substitute_variables__(const RingType& ring, const MapType& idx2poly, NaviType navi) {
00955 
00956   if (navi.isConstant())
00957     return ring.constant(navi.terminalValue());
00958 
00959   //typename MapType::const_reference var(idx2poly[*navi]);
00960   //  PBORI_ASSERT(var.ring() == ring);
00961   return (idx2poly[*navi] * 
00962           substitute_variables__<PolyType>(ring, idx2poly, navi.thenBranch())) + 
00963     substitute_variables__<PolyType>(ring, idx2poly, navi.elseBranch());
00964 } 
00965 
00968 template <class RingType, class MapType, class PolyType>
00969 PolyType
00970 substitute_variables(const RingType& ring,
00971                      const MapType& idx2poly, const PolyType& poly) {
00972 
00973   return substitute_variables__<PolyType>(ring, idx2poly, poly.navigation());
00974 } 
00975 
00976 
00977 END_NAMESPACE_PBORI
00978 
00979 #endif