And now let's code large_bitset.

template 
<
    typename    DomainT = nat64, 
    typename    BitSetT = bits64, 
    ICL_COMPARE Compare = ICL_COMPARE_INSTANCE(std::less, DomainT),
    ICL_INTERVAL(ICL_COMPARE) Interval = ICL_INTERVAL_INSTANCE(ICL_INTERVAL_DEFAULT, DomainT, Compare),
    ICL_ALLOC   Alloc   = std::allocator
> 
class large_bitset
    : boost::equality_comparable < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>
    , boost::less_than_comparable< large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>

    , boost::addable       < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>
    , boost::orable        < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>
    , boost::subtractable  < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>
    , boost::andable       < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>
    , boost::xorable       < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>

    , boost::addable2      < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>, DomainT
    , boost::orable2       < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>, DomainT
    , boost::subtractable2 < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>, DomainT
    , boost::andable2      < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>, DomainT
    , boost::xorable2      < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>, DomainT

    , boost::addable2      < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>, ICL_INTERVAL_TYPE(Interval,DomainT,Compare)
    , boost::orable2       < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>, ICL_INTERVAL_TYPE(Interval,DomainT,Compare)
    , boost::subtractable2 < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>, ICL_INTERVAL_TYPE(Interval,DomainT,Compare)
    , boost::andable2      < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>, ICL_INTERVAL_TYPE(Interval,DomainT,Compare)
    , boost::xorable2      < large_bitset<DomainT,BitSetT,Compare,Interval,Alloc>, ICL_INTERVAL_TYPE(Interval,DomainT,Compare)
      > > > > > > > > > > > > > > > > >
    //^ & - | + ^ & - | + ^ & - | + < == 
    //segment   element   container

The first template parameter DomainT will be instantiated with an integral type that defines the kind of numbers that can be elements of the set. Since we want to go for a large set we use nat64 as default which is a 64 bit unsigned integer ranging from 0 to 2^64-1. As bitset parameter we also choose a 64-bit default. Parameters Combine and Interval are necessary to be passed to dependent type expressions. An allocator can be chosen, if desired.

The nested list of private inheritance contains groups of template instantiations from Boost.Operator, that provides derivable operators from more fundamental once. Implementing the fundamental operators, we get the derivable ones for free. Below is a short overview of what we get using Boost.Operator, where S stands for large_bitset, i for it's interval_type and e for it's domain_type or element_type.

There is a number of associated types

typedef boost::icl::interval_map
    <DomainT, BitSetT, boost::icl::partial_absorber, 
     std::less, boost::icl::inplace_bit_add, boost::icl::inplace_bit_and> interval_bitmap_type;

typedef DomainT                                      domain_type;
typedef DomainT                                      element_type;
typedef BitSetT                                      bitset_type;
typedef typename BitSetT::word_type                  word_type;
typedef typename interval_bitmap_type::interval_type interval_type;
typedef typename interval_bitmap_type::value_type    value_type;

most importantly the implementing interval_bitmap_type that is used for the implementing container.

private:
    interval_bitmap_type _map;

In order to use Boost.Operator we have to implement the fundamental operators as class members. This can be done quite schematically.

public:
    bool     operator ==(const large_bitset& rhs)const { return _map == rhs._map; }
    bool     operator < (const large_bitset& rhs)const { return _map <  rhs._map; }

    large_bitset& operator +=(const large_bitset& rhs) {_map += rhs._map; return *this;}
    large_bitset& operator |=(const large_bitset& rhs) {_map |= rhs._map; return *this;}
    large_bitset& operator -=(const large_bitset& rhs) {_map -= rhs._map; return *this;}
    large_bitset& operator &=(const large_bitset& rhs) {_map &= rhs._map; return *this;}
    large_bitset& operator ^=(const large_bitset& rhs) {_map ^= rhs._map; return *this;}

    large_bitset& operator +=(const element_type& rhs) {return add(interval_type(rhs));      }
    large_bitset& operator |=(const element_type& rhs) {return add(interval_type(rhs));      }
    large_bitset& operator -=(const element_type& rhs) {return subtract(interval_type(rhs)); }
    large_bitset& operator &=(const element_type& rhs) {return intersect(interval_type(rhs));}
    large_bitset& operator ^=(const element_type& rhs) {return flip(interval_type(rhs));     }

    large_bitset& operator +=(const interval_type& rhs){return add(rhs);      }
    large_bitset& operator |=(const interval_type& rhs){return add(rhs);      }
    large_bitset& operator -=(const interval_type& rhs){return subtract(rhs); }
    large_bitset& operator &=(const interval_type& rhs){return intersect(rhs);}
    large_bitset& operator ^=(const interval_type& rhs){return flip(rhs);     }
    

As we can see, the seven most important operators that work on the class type large_bitset can be directly implemented by propagating the operation to the implementing _map of type interval_bitmap_type. For the operators that work on segment and element types, we use member functions add, subtract, intersect and flip. As we will see only a small amount of adaper code is needed to couple those functions with the functionality of the implementing container.

Member functions add, subtract, intersect and flip, that allow to combine intervals to large_bitsets can be uniformly implemented using a private function segment_apply that applies addition, subtraction, intersection or symmetric difference, after having translated the interval's borders into the right bitset positions.

large_bitset& add      (const interval_type& rhs){return segment_apply(&large_bitset::add_,      rhs);}
large_bitset& subtract (const interval_type& rhs){return segment_apply(&large_bitset::subtract_, rhs);}
large_bitset& intersect(const interval_type& rhs){return segment_apply(&large_bitset::intersect_,rhs);}
large_bitset& flip     (const interval_type& rhs){return segment_apply(&large_bitset::flip_,     rhs);}

In the sample programs, that we will present to demonstrate the capabilities of large_bitset we will need a few additional functions specifically output functions in two different flavors.

size_t interval_count()const { return boost::icl::interval_count(_map); }

void show_segments()const
{
    for(typename interval_bitmap_type::const_iterator it_ = _map.begin();
        it_ != _map.end(); ++it_)
    {
        interval_type   itv  = it_->first;
        bitset_type     bits = it_->second;
        std::cout << itv << "->" << bits.as_string("01") << std::endl;
    }
}

void show_matrix(const char off_on[2] = " 1")const
{
    using namespace boost;
    typename interval_bitmap_type::const_iterator iter = _map.begin();
    while(iter != _map.end())
    {
        element_type fst = icl::first(iter->first), lst = icl::last(iter->first);
        for(element_type chunk = fst; chunk <= lst; chunk++)
            std::cout << iter->second.as_string(off_on) << std::endl;
        ++iter;
    }
}

In order to implement operations like the addition of an element say 42 to the large bitset, we need to translate the value to the position of the associated bit representing 42 in the interval container of bitsets. As an example, suppose we use a

large_bitset<nat, mini::bits8> lbs;

that carries small bitsets of 8 bits only. The first four interval of lbs are assumed to be associated with some bitsets. Now we want to add the interval [a,b]==[5,27]. This will result in the following situation:

   [0,1)->   [1,2)->   [2,3)->   [3,4)->   
   [00101100][11001011][11101001][11100000]
+       [111  11111111  11111111  1111]      [5,27] as bitset
         a                           b
         
=> [0,1)->   [1,3)->   [3,4)->
   [00101111][11111111][11110000]

So we have to convert values 5 and 27 into a part that points to the interval and a part that refers to the position within the interval, which is done by a division and a modulo computation. (In order to have a consistent representation of the bitsequences across the containers, within this project, bitsets are denoted with the least significant bit on the left!)

A = a/8 =  5/8 = 0 // refers to interval 
B = b/8 = 27/8 = 3
R = a%8 =  5%8 = 5 // refers to the position in the associated bitset.
S = b%8 = 27%8 = 3

All division and modulo operations needed here are always done using a divisor d that is a power of 2: d = 2^x. Therefore division and modulo can be expressed by bitset operations. The constants needed for those bitset computations are defined here:

private:                                      // Example value
    static const word_type                    //   8-bit case  
        digits  = std::numeric_limits         // --------------------------------------------------------------
                  <word_type>::digits       , //   8           Size of the associated bitsets 
        divisor = digits                    , //   8           Divisor to find intervals for values
        last    = digits-1                  , //   7           Last bit (0 based)
        shift   = log2_<divisor>::value     , //   3           To express the division as bit shift
        w1      = static_cast<word_type>(1) , //               Helps to avoid static_casts for long long
        mask    = divisor - w1              , //   7=11100000  Helps to express the modulo operation as bit_and
        all     = ~static_cast<word_type>(0), // 255=11111111  Helps to express a complete associated bitset
        top     = w1 << (digits-w1)      ;    // 128=00000001  Value of the most significant bit of associated bitsets
                                              //            !> Note: Most signigicant bit on the right.
    

Looking at the example again, we can see that we have to identify the positions of the beginning and ending of the interval [5,27] that is to insert, and then subdivide that range of bitsets into three partitions.

combine interval [5,27] to large_bitset lbs w.r.t. some operation o

   [0,1)->   [1,2)->   [2,3)->   [3,4)->   
   [00101100][11001011][11101001][11100000]
o       [111  11111111  11111111  1111]
         a                           b
subdivide:
   [first!  ][mid_1] . . .[mid_n][   !last]
   [00000111][1...1] . . .[1...1][11110000]

After subdividing, we perform the operation o as follows:

The algorithm, that has been outlined and illustrated by the example, is implemented by the private member function segment_apply. To make the combiner operation a variable in this algorithm, we use a pointer to member function type

typedef void (large_bitset::*segment_combiner)(element_type, element_type, bitset_type);

as first function argument. We will pass member functions combine_ here,

combine_(first_of_interval, end_of_interval, some_bitset);

that take the beginning and ending of an interval and a bitset and combine them to the implementing interval_bitmap_type _map. Here are these functions:

void       add_(DomainT lo, DomainT up, BitSetT bits){_map += value_type(interval_type::right_open(lo,up), bits);}
void  subtract_(DomainT lo, DomainT up, BitSetT bits){_map -= value_type(interval_type::right_open(lo,up), bits);}
void intersect_(DomainT lo, DomainT up, BitSetT bits){_map &= value_type(interval_type::right_open(lo,up), bits);}
void      flip_(DomainT lo, DomainT up, BitSetT bits){_map ^= value_type(interval_type::right_open(lo,up), bits);}

Finally we can code function segment_apply, that does the partitioning and subsequent combining:

large_bitset& segment_apply(segment_combiner combine, const interval_type& operand)
{
    using namespace boost;
    if(icl::is_empty(operand))
        return *this;
                                                        // same as
    element_type   base = icl::first(operand) >> shift, // icl::first(operand) / divisor
                   ceil = icl::last (operand) >> shift; // icl::last (operand) / divisor
    word_type base_rest = icl::first(operand) &  mask , // icl::first(operand) % divisor
              ceil_rest = icl::last (operand) &  mask ; // icl::last (operand) % divisor  

    if(base == ceil) // [first, last] are within one bitset (chunk)
        (this->*combine)(base, base+1, bitset_type(  to_upper_from(base_rest)
                                                   & from_lower_to(ceil_rest)));
    else // [first, last] spread over more than one bitset (chunk)
    {
        element_type mid_low = base_rest == 0   ? base   : base+1, // first element of mid part 
                     mid_up  = ceil_rest == all ? ceil+1 : ceil  ; // last  element of mid part

        if(base_rest > 0)    // Bitset of base interval has to be filled from base_rest to last
            (this->*combine)(base, base+1, bitset_type(to_upper_from(base_rest)));
        if(ceil_rest < all)  // Bitset of ceil interval has to be filled from first to ceil_rest
            (this->*combine)(ceil, ceil+1, bitset_type(from_lower_to(ceil_rest)));
        if(mid_low < mid_up) // For the middle part all bits have to set.
            (this->*combine)(mid_low, mid_up, bitset_type(all));
    }
    return *this;
}

The functions that help filling bitsets to and from a given bit are implemented here:

static word_type from_lower_to(word_type bit){return bit==last ? all : (w1<<(bit+w1))-w1;}
static word_type to_upper_from(word_type bit){return bit==last ? top : ~((w1<<bit)-w1); }

This completes the implementation of class template large_bitset. Using only a small amount of mostly schematic code, we have been able to provide a pretty powerful, self compressing and generally usable set type for all integral domain types.