main page namespaces classes files
ELL home

src/ell/Basin.cc

Go to the documentation of this file.

#include "ell/Basin.hh"
#include "ell/util/PriorityQueue.hh"
#include <biu/assertbiu.hh>
#include <list>


namespace ell
{

    Basin::QueueVal::QueueVal()
     :  minNeigh(), minNeighE(0.0), minNeighGWE(UINT_MAX)
    {
    }
    
    Basin::QueueVal::~QueueVal()
    {
    }


    Basin::Basin(   const double maxE_, 
                    const double deltaE_, 
                    const size_t maxToStore_)
      : maxE(maxE_),
        deltaE(deltaE_),
        maxToStore(maxToStore_)
    {}
    
    double
    Basin::flood(   const ell::State* const localMin, 
                    StateCollector* scBasin,
                    StateCollector* scSurface)
    {
        return Basin::basin( localMin, scBasin, 
                                maxE, deltaE, maxToStore, scSurface);
    }

    double 
    Basin::basin(   const State* const localMin, 
            StateCollector* scBasin,
            const double maxE_, 
            const double deltaE, 
            const size_t maxToStore, 
            StateCollector* scSurface
            ) 
    {
        assertbiu( localMin != NULL, "No local Min given! (==NULL)"); 
        assertbiu( scBasin != NULL, "No StateCollector given! (scBasin==NULL)"); 
        assertbiu( localMin->getEnergy() < maxE_, "E(localMin) >= maxE");
        
        using namespace ell::util;
        
          // the priority queue to store neighbors to handle in
        PriorityQueue<QueueVal> pq;
          // init using the given local minimum
        PriorityQueue<QueueVal>::iterator it = pq.insert(localMin);
        it->second.minNeigh.clear(); // has no minimal neighbor (is min)
        it->second.minNeighE = it->first.E;
        it->second.minNeighGWE = 0;  
        
        double maxE = maxE_;
        
        typedef std::list< PriorityQueue<QueueVal>::QueueKey  > NeighList;
        
        while (!pq.empty()) {
              // get top element
            PriorityQueue<QueueVal>::const_iterator top = pq.top();
              // uncompress top element into new State
            State* s = localMin->uncompress(top->first.s, NULL);
            CSequence minNeigh; // neighbor with minimal energy
            double minNeighE = (double)INT_MAX; // energy of the minimal neighbor
             
            
              // get neighbor list
            State::NeighborListPtr neighbors = s->getNeighborList();
            State::NeighborList::Iterator it = neighbors->begin();
                // contains all neighbors with E(top) < E < maxE
            NeighList toStore;
              // iterate through all neighbors
            while (neighbors->end() != it) 
            {
                  // get neighbor data
                double nextE = it->getEnergy();
                CSequence next = it->compress();
                  // increase iterator for next iteration
                ++it;

                  // check for lowest neighbor
                if (    (nextE < minNeighE) ||
                        (nextE == minNeighE && next < minNeigh ) ) 
                {
                    minNeigh = next;
                    minNeighE = nextE;
                }
                  // neighbor rejection
                if (    nextE < maxE    // energy does not exceed maxE
                        && ( nextE > top->first.E  // energy higher than top or
                        || (nextE == top->first.E && next > top->first.s) 
                                // greater in sequence order 
                        ) ) 
                {
                    // store neighbor with higher energy that has to be added to pq
                    toStore.push_front(NeighList::value_type(nextE, next));
                }
                    
            }
              // handle top depending on its basin membership           
            if (    (minNeighE == top->second.minNeighE 
                    && minNeigh == top->second.minNeigh)   // belongs to the basin
                    || top->second.minNeigh.empty() ) // top is the local minimum
            {
                  // store because top belongs to the basin 
                scBasin->add(*s);
                
                PriorityQueue<QueueVal>::iterator toFill;
                  // sort list with increasing energy
                toStore.sort();
                for (   NeighList::const_iterator n = toStore.begin(); 
                        n != toStore.end() && n->E < maxE; 
                        n++) 
                {
                    toFill = pq.find(*n);
                    if (toFill == pq.end()) {
                          // insert empty value object
                        toFill = pq.insert(*n);
                          // feed data
                        toFill->second.minNeigh = top->first.s;
                        toFill->second.minNeighE = top->first.E;
                        toFill->second.minNeighGWE = top->second.minNeighGWE;
                          // check queue size and reduce if necessary
                        if ( pq.size() >= maxToStore ) { // queue to big
                              // get current maximal E in queue
                            double curMaxE = pq.getMaxE();
                              // reduce maximally allowed energy below cur. max
                            do { maxE -= deltaE; } while (maxE > curMaxE);
                              // erase all elements >= maxE from the queue;
                            pq.reduceMaxE(maxE);
                        }
                    } else  // update of entry if this one is smaller
                        if (    toFill->second.minNeighE > top->first.E ||
                                ( toFill->second.minNeighE == top->first.E &&
                                  toFill->second.minNeigh > top->first.s) )
                    {
                          // feed data
                        toFill->second.minNeigh = top->first.s;
                        toFill->second.minNeighE = top->first.E;
                        toFill->second.minNeighGWE = top->second.minNeighGWE;
                    }
                }
                
            } else { // belongs NOT to the basin but is surface
                if (scSurface != NULL) {
                      // store as processed state 
                    scSurface->add(*s);
                }
            } 

              // clear memory
            delete s;
              // remove top element
            pq.pop();           
        }
         
        return maxE;
    }


}