LimitedCombination.hpp

Go to the documentation of this file.

#ifndef __NUMERIC_LIMITED_COMBINATION_HPP__
#define __NUMERIC_LIMITED_COMBINATION_HPP__

#include <map>
#include <vector>
#include <numeric/Combinatorics.hpp>

namespace numeric {

template<typename AtomType>
class LimitedCombination
{
    typedef uint32_t CoreType;
    typedef std::map<AtomType, size_t> AtomType2CountMap;
    typedef std::vector<AtomType> AtomTypeList;
    typedef std::vector<CoreType> CombinationDescriptor;

    // Model with the number of available models of this type
    AtomType2CountMap mAtomTypeAvailablilityMap;
    AtomTypeList mAtomTypeList;

    std::vector<CoreType> mLimits;
    std::vector<CoreType> mCurrentCombination;

    size_t mSize;
    numeric::Mode mMode;

    void prepare()
    {
        CoreType currentType = 0;
        typename AtomType2CountMap::const_iterator cit = mAtomTypeAvailablilityMap.begin();
        for(; cit != mAtomTypeAvailablilityMap.end(); ++cit, ++currentType)
        {
            mAtomTypeList.push_back(cit->first);
            size_t numberPerType = cit->second;
            mLimits.at(currentType) = numberPerType;
        }

        // Making sure this behaves like the other combinatorics loops
        if(!next())
        {
            throw std::runtime_error("numeric::LimitedCombination::prepare: preparation failed "
                    "check the given parameters");
        }
    }

    std::vector<AtomType> mapToAtomTypes(const std::vector<CoreType>& combination) const
    {
        std::vector<AtomType> atomTypeList;
        for(size_t coreType = 0; coreType < combination.size(); ++coreType)
        {
            size_t count = combination[coreType];
            for(size_t i = 0; i < count; ++i)
            {
                atomTypeList.push_back( mAtomTypeList[coreType] );
            }
        }
        return atomTypeList;
    }

    bool increment(std::vector<uint32_t>& combination, const std::vector<uint32_t>& limits) const
    {
        for(size_t pos = 0; pos < combination.size(); ++pos)
        {
            if( combination[pos] < limits[pos] )
            {
                ++combination[pos];
                if(pos == 0)
                {
                    return true;
                }

                // reset lower values after overflow
                for(int lpos = pos-1; lpos >= 0; --lpos)
                {
                    combination[lpos] = 0;
                }
                return true;
            } else if(pos == combination.size() - 1)
            {
                return false;
            }
        }
        throw std::runtime_error("numeric::LimitedCombination::increment internal error -- "
                " please check consistency of your input");
    }

public:
    LimitedCombination(const AtomType2CountMap& countMap, size_t size, numeric::Mode mode)
        : mAtomTypeAvailablilityMap(countMap)
        , mLimits(countMap.size(),0)
        , mCurrentCombination(countMap.size(),0)
        , mSize(size)
        , mMode(mode)
    {
        size_t totalCount = totalNumberOfAtoms(mAtomTypeAvailablilityMap);
        if(countMap.empty() || totalCount == 0)
        {
            throw std::invalid_argument("numeric::LimitedCombination no atoms to generate combination from -- check for empty map");
        }
        if(mSize > totalCount)
        {
            mSize = totalCount;
        }

        prepare();
    }

    static size_t totalNumberOfAtoms(const AtomType2CountMap& countMap)
    {
        typename AtomType2CountMap::const_iterator cit = countMap.begin();
        size_t count = 0;
        for(; cit != countMap.end(); ++cit)
        {
            count += cit->second;
        }
        return count;
    }

    std::vector<AtomType> current() const
    {
        std::vector<AtomType> atomTypeList = mapToAtomTypes(mCurrentCombination);
        std::sort(atomTypeList.begin(), atomTypeList.end());
        return atomTypeList;
    }

    size_t getCombinationSize(const std::vector<CoreType>& combination) const
    {
        size_t combinationSize = 0;
        std::vector<CoreType>::const_iterator cit = combination.begin();
        for(; cit != combination.end(); ++cit)
        {
            combinationSize += *cit;
        }
        return combinationSize;
    }

    bool next()
    {
        while(true)
        {
            bool hasNext = increment(mCurrentCombination, mLimits);
            if(hasNext)
            {
                // check for size constraint
                size_t combinationSize = getCombinationSize(mCurrentCombination);

                switch(mMode)
                {
                    case EXACT:
                        if(combinationSize != mSize)
                        {
                            continue;
                        }
                        break;
                    case MAX:
                        if(combinationSize > mSize)
                        {
                            continue;
                        }
                        break;
                    case MIN:
                        if(combinationSize < mSize)
                        {
                            continue;
                        }
                        break;
                }
                return true;
            }

            return false;
        }

        return false;
    }

};

} // numeric
#endif // __NUMERIC_LIMITED_COMBINATION_HPP__